High level central plant optimization

ABSTRACT

A controller for equipment that operate to provide heating or cooling to a building or campus includes a processing circuit configured to obtain utility rate data indicating a price of resources consumed by the equipment to serve energy loads of the building or campus, obtain an objective function that expresses a total monetary cost of operating the equipment over an optimization period as a function of the utility rate data and an amount of the resources consumed by the equipment, determine a relationship between resource consumption and load production of the equipment, optimize the objective function over the optimization subject to a constraint based on the relationship between the resource consumption and the load production of the equipment to determine a distribution of the load production across the equipment, and operate the equipment to achieve the distribution.

CROSS-REFERENCE TO RELATED PATENT APPLICATIONS

This application is a continuation of U.S. patent application Ser. No.16/214,984 filed Dec. 10, 2018, which is a continuation of U.S. patentapplication Ser. No. 14/634,609 filed Feb. 27, 2015, and granted as U.S.Pat. No. 10,175,681 on Jan. 8, 2019, which claims the benefit of andpriority to U.S. Provisional Patent Application No. 61/987,361 filed May1, 2014. The entire disclosures of each of these patent applications andpatent are incorporated by reference herein.

BACKGROUND

The present disclosure relates generally to the operation of a centralplant for serving building thermal energy loads. The present disclosurerelates more particularly to systems and methods for distributingbuilding thermal energy loads across a plurality of subplants configuredto serve the building thermal energy loads.

A central plant may include various types of equipment configured toserve the thermal energy loads of a building or campus (i.e., a systemof buildings). For example, a central plant may include heaters,chillers, heat recovery chillers, cooling towers, or other types ofequipment configured to provide heating or cooling for the building. Acentral plant may consume resources from a utility (e.g., electricity,water, natural gas, etc.) to heat or cool a working fluid (e.g., water,glycol, etc.) that is circulated to the building or stored for later useto provide heating or cooling for the building. Fluid conduits typicallydeliver the heated or chilled fluid to air handlers located on therooftop of the building or to individual floors or zones of thebuilding. The air handlers push air past heat exchangers (e.g., heatingcoils or cooling coils) through which the working fluid flows to provideheating or cooling to the air. The working fluid then returns to thecentral plant to receive further heating or cooling and the cyclecontinues.

High efficiency equipment can help reduce the amount of energy consumedby a central plant; however, the effectiveness of such equipment ishighly dependent on the control technology that is used to distributethe load across the multiple subplants. For example, it may be more costefficient to run heat pump chillers instead of conventional chillers anda water heater when energy prices are high. It is difficult andchallenging to determine when and to what extent each of the multiplesubplants should be used to minimize energy cost. If electrical demandcharges are considered, the optimization is even more complicated.

Thermal energy storage can be used to store energy for later use. Whencoupled with real-time pricing for electricity and demand charges,thermal energy storage provides a degree of flexibility that can be usedto greatly decrease energy costs by shifting production to low costtimes or when other electrical loads are lower so that a new peak demandis not set. It is difficult and challenging to integrate thermal energystorage with a central plant having multiple subplants and to optimizethe use of thermal energy storage in conjunction with the multiplesubplants to minimize energy cost.

SUMMARY

One implementation of the present disclosure is an optimization systemfor a central plant configured to serve building energy loads. Theoptimization system includes a processing circuit configured to receiveload prediction data indicating building energy loads for a plurality oftime steps in an optimization period and utility rate data indicating aprice of one or more resources consumed by equipment of the centralplant to serve the building energy loads at each of the plurality oftime steps. The optimization system further includes a high leveloptimization module configured to generate an objective function thatexpresses a total monetary cost of operating the central plant over theoptimization period as a function of the utility rate data and an amountof the one or more resources consumed by the central plant equipment ateach of the plurality of time steps. The high level optimization moduleis configured to optimize the objective function over the optimizationperiod subject to load equality constraints and capacity constraints onthe central plant equipment to determine an optimal distribution of thebuilding energy loads over multiple groups of the central plantequipment at each of the plurality of time steps.

In some embodiments, the high level optimization module uses linearprogramming to generate and optimize the objective function. Theobjective function may include a cost vector having cost variablesrepresenting a monetary cost associated with each of the one or moreresources consumed by the central plant equipment to serve the buildingenergy loads at each of the plurality of time steps. The objectivefunction may further include a decision matrix comprising load variablesrepresenting an energy load for each of the multiple groups of thecentral plant equipment at each of the plurality of time steps. The highlevel optimization module may be configured to determine optimal valuesfor the load variables in the decision matrix.

In some embodiments, the central plant includes a plurality ofsubplants. Each of the multiple groups of the central plant equipmentmay correspond to one of the plurality of subplants. In someembodiments, the plurality of subplants include at least one of a hotthermal energy storage subplant and a cold thermal energy storagesubplant. The thermal energy storage subplants may be configured tostore thermal energy generated in one of the plurality of time steps foruse in another of the plurality of time steps.

In some embodiments, the high level optimization module is configured togenerate a subplant curve for each of the plurality of subplants. Eachsubplant curve may indicate a relationship between resource consumptionand load production for one of the plurality of subplants. The highlevel optimization module may use the subplant curves to formulatesubplant curve constraints and may optimize the objective functionsubject to the subplant curve constraints. In some embodiments,generating the subplant curve includes at least one of converting anonlinear subplant curve into a linear subplant curve including one ormore piecewise linear segments and converting a non-convex subplantcurve into a convex subplant curve. In some embodiments, generating thesubplant curve includes receiving an initial subplant curve based onmanufacturer data for the group of equipment corresponding to thesubplant and updating the initial subplant curve using experimental datafrom the central plant.

Another implementation of the present disclosure is a cascadedoptimization system for a central plant configured to serve buildingenergy loads. The cascaded optimization system includes a central plantcontroller configured to use dynamic programming to split anoptimization problem for the central plant into a high leveloptimization and a low level optimization, a high level optimizationmodule configured to perform the high level optimization, and a lowlevel optimization module configured to perform the low leveloptimization. The high level optimization includes determining anoptimal distribution of building energy loads over multiple groups ofcentral plant equipment. The low level optimization includes determiningoptimal operating statuses for individual devices within each of themultiple groups of the central plant equipment.

In some embodiments, the optimal distribution of the building energyloads determined by the high level optimization module optimizes amonetary cost of operating the central plant over an optimizationperiod. The optimal operating statuses determined by the low leveloptimization module may optimize an amount of energy consumed by each ofthe multiple groups of the central plant equipment to achieve theoptimal distribution of the building energy loads determined by the highlevel optimization module.

In some embodiments, the low level optimization module is configured togenerate a subplant curve for each of the groups of central plantequipment. Each subplant curve may indicate a relationship betweenresource consumption and load production for one of the groups ofcentral plant equipment. The high level optimization module may beconfigured to use the subplant curves to formulate subplant curveconstraints and to determine the optimal distribution of the buildingenergy loads subject to the subplant curve constraints.

Another implementation of the present disclosure is a method foroptimizing cost in a central plant configured to serve building energyloads. The method includes receiving, at a processing circuit of acentral plant optimization system, load prediction data indicatingbuilding energy loads for a plurality of time steps in an optimizationperiod and utility rate data indicating a price of one or more resourcesconsumed by equipment of the central plant to serve the building energyloads at each of the plurality of time steps. The method furtherincludes generating, by a high level optimization module of the centralplant optimization system, an objective function that expresses a totalmonetary cost of operating the central plant over the optimizationperiod as a function of the utility rate data and an amount of the oneor more resources consumed by the central plant equipment at each of theplurality of time steps. The method further includes optimizing, by thehigh level optimization module, the objective function over theoptimization period subject to load equality constraints and capacityconstraints on the central plant equipment to determine an optimaldistribution of the building energy loads over multiple groups of thecentral plant equipment at each of the plurality of time steps.

In some embodiments, the high level optimization module uses linearprogramming to generate and optimize the objective function. Theobjective function may include a cost vector having cost variablesrepresenting a monetary cost associated with each of the one or moreresources consumed by the central plant equipment to serve the buildingenergy loads at each of the plurality of time steps. The objectivefunction may further include a decision matrix comprising load variablesrepresenting an energy load for each of the multiple groups of thecentral plant equipment at each of the plurality of time steps. The highlevel optimization module may be configured to determine optimal valuesfor the load variables in the decision matrix.

In some embodiments, the method includes using the building energy loadsand capacity limits for the central plant equipment to generate the loadequality constraints and the capacity constraints. The load equalityconstraints may ensure that the optimal distribution satisfies thebuilding energy loads at each of the plurality of time steps. Thecapacity constraints may ensure that the multiple groups of centralplant equipment are operated within the capacity limits at each of theplurality of time steps.

In some embodiments, the central plant includes a plurality ofsubplants. Each of the multiple groups of the central plant equipmentmay correspond to one of the plurality of subplants. In someembodiments, the method includes generating a subplant curve for each ofthe plurality of subplants. Each subplant curve may indicate arelationship between resource consumption and load production for one ofthe plurality of subplants. The method may further include using thesubplant curves to formulate subplant curve constraints and optimizingthe objective function subject to the subplant curve constraints. Insome embodiments, generating the subplant curve includes at least one ofconverting a nonlinear subplant curve into a linear subplant curveincluding one or more piecewise linear segments and converting anon-convex subplant curve into a convex subplant curve.

In some embodiments, the method further includes using dynamicprogramming to split the optimization process into a high leveloptimization and a low level optimization. The high level optimizationmay include determining the optimal distribution of the building energyloads over the multiple groups of the central plant equipment. In someembodiments, the optimal distribution of the building energy loadsoptimizes the monetary cost of operating the central plant over theoptimization period. The low level optimization may include determiningoptimal operating statuses for individual devices within each of themultiple groups of the central plant equipment. In some embodiments, theoptimal operating statuses optimize an amount of energy consumed by eachof the multiple groups of the central plant equipment to achieve theoptimal distribution of the building energy loads.

Those skilled in the art will appreciate that the summary isillustrative only and is not intended to be in any way limiting. Otheraspects, inventive features, and advantages of the devices and/orprocesses described herein, as defined solely by the claims, will becomeapparent in the detailed description set forth herein and taken inconjunction with the accompanying drawings.

BRIEF DESCRIPTION OF THE DRAWINGS

FIG. 1 is a schematic diagram of a central plant having a plurality ofsubplants including a heater subplant, heat recovery chiller subplant, achiller subplant, a hot thermal energy storage subplant, and a coldthermal energy storage subplant, according to an exemplary embodiment.

FIG. 2 is a block diagram illustrating a central plant system includinga central plant controller that may be used to control the central plantof FIG. 1 , according to an exemplary embodiment.

FIG. 3 is block diagram illustrating a portion of central plant systemof FIG. 2 in greater detail, showing a load/rate prediction module, ahigh level optimization module, a low level optimization module, abuilding automation system, and central plant equipment, according to anexemplary embodiment.

FIG. 4 , a block diagram illustrating the high level optimization moduleof FIG. 3 in greater detail, according to an exemplary embodiment.

FIGS. 5A-5B are subplant curves illustrating a relationship between theresource consumption of a subplant and the subplant load and which maybe used by the high level optimization module of FIG. 4 to optimize theperformance of the central plant of FIG. 1 , according to an exemplaryembodiment.

FIG. 6 is a non-convex and nonlinear subplant curve that may begenerated from experimental data or by combining equipment curves forindividual devices of the central plant, according to an exemplaryembodiment.

FIG. 7 is a linearized subplant curve that may be generated from thesubplant curve of FIG. 6 by converting the non-convex and nonlinearsubplant curve into piecewise linear segments, according to an exemplaryembodiment.

FIG. 8 is a graph illustrating a set of subplant curves that may begenerated by the high level optimization module of FIG. 3 based onexperimental data from a low level optimization module for multipledifferent environmental conditions, according to an exemplaryembodiment.

FIG. 9 is a block diagram of a planning system that incorporates thehigh level optimization module of FIG. 3 , according to an exemplaryembodiment.

FIG. 10 is a drawing illustrating the operation of the planning systemof FIG. 9 , according to an exemplary embodiment.

FIG. 11 is a flowchart of a process for optimizing cost in a centralplant that may be performed by the central plant controller of FIG. 2 orthe planning system of FIG. 9 , according to an exemplary embodiment.

DETAILED DESCRIPTION

Referring generally to the FIGURES, systems and methods for optimizing acentral plant are shown, according to an exemplary embodiment. A centralplant may include may include various types of equipment configured toserve the thermal energy loads of a building or campus (i.e., a systemof buildings). For example, a central plant may include heaters,chillers, heat recovery chillers, cooling towers, or other types ofequipment configured to provide heating or cooling for the building orcampus. The central plant equipment may be divided into various groupsconfigured to perform a particular function. Such groups of centralplant equipment are referred to herein as subplants. For example, acentral plant may include a heater subplant, a chiller subplant, a heatrecovery chiller subplant, a cold thermal energy storage subplant, a hotthermal energy storage subplant, etc. The subplants may consumeresources from one or more utilities (e.g., water, electricity, naturalgas, etc.) to serve the energy loads of the building or campus.Optimizing the central plant may include operating the various subplantsin such a way that results in a minimum monetary cost to serve thebuilding energy loads.

In some embodiments, the central plant optimization is a cascadedoptimization process including a high level optimization and a low leveloptimization. The high level optimization may determine an optimaldistribution of energy loads across the various subplants. For example,the high level optimization may determine a thermal energy load to beproduced by each of the subplants at each time element in anoptimization period. In some embodiments, the high level optimizationincludes optimizing a high level cost function that expresses themonetary cost of operating the subplants as a function of the resourcesconsumed by the subplants at each time element of the optimizationperiod. The low level optimization may use the optimal load distributiondetermined by the high level optimization to determine optimal operatingstatuses for individual devices within each subplant. Optimal operatingstatuses may include, for example, on/off states and/or operatingsetpoints for individual devices of each subplant. The low leveloptimization may include optimizing a low level cost function thatexpresses the energy consumption of a subplant as a function of theon/off states and/or operating setpoints for the individual devices ofthe subplant.

The present disclosure focuses on the high level optimization anddescribes systems and methods for performing the high leveloptimization. A high level optimization module may perform the highlevel optimization. In various embodiments, the high level optimizationmodule may be a component of a central plant controller configured forreal-time control of a physical plant or a component of a planning toolconfigured to optimize a simulated plant (e.g., for planning or designpurposes).

In some embodiments, the high level optimization module uses a linearprogramming framework to perform the high level optimization.Advantageously, linear programming can efficiently handle complexoptimization scenarios and can optimize over a relatively longoptimization period (e.g., days, weeks, years, etc.) in a relativelyshort timeframe (e.g., seconds, milliseconds, etc.). In otherembodiments, the high level optimization module may use any of a varietyof other optimization frameworks (e.g., quadratic programming,linear-fractional programming, nonlinear programming, combinatorialalgorithms, etc.).

An objective function defining the high level optimization problem canbe expressed in the linear programming framework as:

${{\underset{x}{\arg\min}c^{T}x};{{{subject}{to}{Ax}} \leq b}},{{Hx} = g}$where c is a cost vector, x is a decision matrix, A and b are a matrixand vector (respectively) which describe inequality constraints on thevariables in the decision matrix x, and H and g are a matrix and vector(respectively) which describe equality constraints on the variables inthe decision matrix x. The variables in the decision matrix x mayinclude the subplant loads assigned to the various subplants and/or anamount of resource consumption by the subplants at each time element inthe optimization period. The high level optimization module may definethe cost vector c and the optimization constraints (e.g., the matrices Aand H and the vectors b and g) and solve the optimization problem todetermine optimal subplant load values for the variables in the decisionmatrix x.

The high level optimization module may receive, as an input, predictedor planned energy loads for the building or campus for each of the timeelements in the optimization period. The high level optimization modulemay use the predicted or planned loads to formulate the constraints onthe high level optimization problem (e.g., to define the matrices A andH and the vectors b and g). The high level optimization module may alsoreceive utility rates (e.g., energy prices, water prices, demandcharges, etc.) defining the cost of each resource consumed by thecentral plant to serve the energy loads. The utility rates may betime-variable rates (e.g., defining a different rates at differenttimes) and may include demand charges for various time periods. The highlevel optimization module may use the utility rates to define the costvector c.

The high level optimization module may receive or generate subplantcurves for each of the subplants. A subplant curve defines the resourceconsumption of a subplant as a function of the load produced by thesubplant. The subplant curves may be generated by a low leveloptimization module or by the high level optimization module based onoperating data points received from the low level optimization module.The high level optimization module may use the subplant curves toconstrain the resource consumption of each subplant to a value along thecorresponding subplant curve (e.g., based on the load produced by thesubplant). For example, the high level optimization module may use thesubplant curves to define the optimization constraints (e.g., thematrices A and H and the vectors b and g) on the high level optimizationproblem.

In some embodiments, the high level optimization module is configured toincorporate a demand charge into the high level optimization process.The demand charge is an additional charge imposed by some utilityproviders based on the maximum rate of resource consumption during anapplicable demand charge period. For example, an electric demand chargemay be provided as a cost c_(demand) per unit power and may bemultiplied by the peak electricity usage max(P_(elec,k)) during a demandcharge period to determine the demand charge. Conventional systems havebeen unable to incorporate a demand charge into a linear optimizationframework due to the nonlinear max( ) function used to calculate thedemand charge.

Advantageously, the high level optimization module of the presentdisclosure may be configured to incorporate the demand charge into thelinear optimization framework by modifying the decision matrix x, thecost vector c, and/or the A matrix and the b vector which describe theinequality constraints. For example, the high level optimization modulemay modify the decision matrix x by adding a new decision variablex_(peak) representing the peak power consumption within the optimizationperiod. The high level optimization module may modify the cost vector cwith the demand charge rate c_(demand) such that the demand charge ratec_(demand) is multiplied by the peak power consumption x_(peak). Thehigh level optimization module may generate and/or impose constraints toensure that the peak power consumption x_(peak) is greater than or equalto the electric demand for each time step in the demand charge periodand greater than or equal to its previous value during the demand chargeperiod.

In some embodiments, the high level optimization module is configured toincorporate a load change penalty into the high level optimizationprocess. The load change penalty may represent an increased cost (e.g.,equipment degradation, etc.) resulting from a rapid change in the loadassigned to a subplant. The high level optimization module mayincorporate the load change penalty by modifying the decision matrix x,the cost vector c, and/or the optimization constraints. For example, thehigh level optimization module may modify the decision matrix x byadding load change variables δ for each subplant. The load changevariables may represent the change in subplant load for each subplantfrom one time element to the next. The high level optimization modulemay modify the cost vector c to add a cost associated with changing thesubplant loads. In some embodiments, the high level optimization moduleadds constraints that constrain the load change variables δ to thecorresponding change in the subplant load. These and other enhancementsto the high level optimization process may be incorporated into thelinear optimization framework, as described in greater detail below.

Referring now to FIG. 1 , a diagram of a central plant 10 is shown,according to an exemplary embodiment. Central plant 10 is shown toinclude a plurality of subplants including a heater subplant 12, a heatrecovery chiller subplant 14, a chiller subplant 16, a cooling towersubplant 18, a hot thermal energy storage (TES) subplant 20, and a coldthermal energy storage (TES) subplant 22. Subplants 12-22 consumeresources (e.g., water, natural gas, electricity, etc.) from utilitiesto serve the thermal energy loads (e.g., hot water, cold water, heating,cooling, etc.) of a building or campus. For example, heater subplant 12may be configured to heat water in a hot water loop 24 that circulatesthe hot water between central plant 10 and a building (not shown).Chiller subplant 16 may be configured to chill water in a cold waterloop 26 that circulates the cold water between central plant 10 and thebuilding. Heat recovery chiller subplant 14 may be configured totransfer heat from cold water loop 26 to hot water loop 24 to provideadditional heating for the hot water and additional cooling for the coldwater. Condenser water loop 28 may absorb heat from the cold water inchiller subplant 16 and reject the absorbed heat in cooling towersubplant 18 or transfer the absorbed heat to hot water loop 24. Hot TESsubplant 20 and cold TES subplant 22 store hot and cold thermal energy,respectively, for subsequent use.

Hot water loop 24 and cold water loop 26 may deliver the heated and/orchilled water to air handlers located on the rooftop of a building or toindividual floors or zones of the building. The air handlers push airpast heat exchangers (e.g., heating coils or cooling coils) throughwhich the water flows to provide heating or cooling for the air. Theheated or cooled air may be delivered to individual zones of thebuilding to serve the thermal energy loads of the building. The waterthen returns to central plant 10 to receive further heating or coolingin subsystems 12-22.

Although central plant 10 is shown and described as heating and coolingwater for circulation to a building, it is understood that any othertype of working fluid (e.g., glycol, CO2, etc.) may be used in place ofor in addition to water to serve the thermal energy loads. In otherembodiments, central plant 10 may provide heating and/or coolingdirectly to the building or campus without requiring an intermediateheat transfer fluid. Central plant 10 may be physically separate from abuilding served by subplants 12-22 or physically integrated with thebuilding (e.g., located within the building).

Each of subplants 12-22 may include a variety of equipment configured tofacilitate the functions of the subplant. For example, heater subplant12 is shown to include a plurality of heating elements 30 (e.g.,boilers, electric heaters, etc.) configured to add heat to the hot waterin hot water loop 24. Heater subplant 12 is also shown to includeseveral pumps 32 and 34 configured to circulate the hot water in hotwater loop 24 and to control the flow rate of the hot water throughindividual heating elements 30. Heat recovery chiller subplant 14 isshown to include a plurality of heat recovery heat exchangers 36 (e.g.,refrigeration circuits) configured to transfer heat from cold water loop26 to hot water loop 24. Heat recovery chiller subplant 14 is also shownto include several pumps 38 and 40 configured to circulate the hot waterand/or cold water through heat recovery heat exchangers 36 and tocontrol the flow rate of the water through individual heat recovery heatexchangers 36.

Chiller subplant 16 is shown to include a plurality of chillers 42configured to remove heat from the cold water in cold water loop 26.Chiller subplant 16 is also shown to include several pumps 44 and 46configured to circulate the cold water in cold water loop 26 and tocontrol the flow rate of the cold water through individual chillers 42.Cooling tower subplant 18 is shown to include a plurality of coolingtowers 48 configured to remove heat from the condenser water incondenser water loop 28. Cooling tower subplant 18 is also shown toinclude several pumps 50 configured to circulate the condenser water incondenser water loop 28 and to control the flow rate of the condenserwater through individual cooling towers 48.

Hot TES subplant 20 is shown to include a hot TES tank 52 configured tostore the hot water for later use. Hot TES subplant 20 may also includeone or more pumps or valves configured to control the flow rate of thehot water into or out of hot TES tank 52. Cold TES subplant 22 is shownto include cold TES tanks 54 configured to store the cold water forlater use. Cold TES subplant 22 may also include one or more pumps orvalves configured to control the flow rate of the cold water into or outof cold TES tanks 54. In some embodiments, one or more of the pumps incentral plant 10 (e.g., pumps 32, 34, 38, 40, 44, 46, and/or 50) orpipelines in central plant 10 includes an isolation valve associatedtherewith. In various embodiments, isolation valves may be integratedwith the pumps or positioned upstream or downstream of the pumps tocontrol the fluid flows in central plant 10. In other embodiments, more,fewer, or different types of devices may be included in central plant10.

Referring now to FIG. 2 , a block diagram illustrating a central plantsystem 100 is shown, according to an exemplary embodiment. System 100 isshown to include a central plant controller 102, a building automationsystem 108, and a plurality of subplants 12-22. Subplants 12-22 may bethe same as previously described with reference to FIG. 1 . For example,subplants 12-22 are shown to include a heater subplant 12, a heatrecovery chiller subplant 14, a chiller subplant 16, a hot TES subplant20, and a cold TES subplant 22.

Each of subplants 12-22 is shown to include equipment 60 that can becontrolled by central plant controller 102 and/or building automationsystem 108 to optimize the performance of central plant 10. Equipment 60may include, for example, heating devices 30, chillers 42, heat recoveryheat exchangers 36, cooling towers 48, thermal energy storage devices52, 54, pumps 32, 44, 50, valves 34, 38, 46, and/or other devices ofsubplants 12-22. Individual devices of equipment 60 can be turned on oroff to adjust the thermal energy load served by each of subplants 12-22.In some embodiments, individual devices of equipment 60 can be operatedat variable capacities (e.g., operating a chiller at 10% capacity or 60%capacity) according to an operating setpoint received from central plantcontroller 102.

In some embodiments, one or more of subplants 12-22 includes a subplantlevel controller configured to control the equipment 60 of thecorresponding subplant. For example, central plant controller 102 maydetermine an on/off configuration and global operating setpoints forequipment 60. In response to the on/off configuration and receivedglobal operating setpoints, the subplant controllers may turn individualdevices of equipment 60 on or off, and implement specific operatingsetpoints (e.g., damper position, vane position, fan speed, pump speed,etc.) to reach or maintain the global operating setpoints.

Building automation system (BAS) 108 may be configured to monitorconditions within a controlled building or building zone. For example,BAS 108 may receive input from various sensors (e.g., temperaturesensors, humidity sensors, airflow sensors, voltage sensors, etc.)distributed throughout the building and may report building conditionsto central plant controller 102. Building conditions may include, forexample, a temperature of the building or a zone of the building, apower consumption (e.g., electric load) of the building, a state of oneor more actuators configured to affect a controlled state within thebuilding, or other types of information relating to the controlledbuilding. BAS 108 may operate subplants 12-22 to affect the monitoredconditions within the building and to serve the thermal energy loads ofthe building.

BAS 108 may receive control signals from central plant controller 102specifying on/off states and/or setpoints for equipment 60. BAS 108 maycontrol equipment 60 (e.g., via actuators, power relays, etc.) inaccordance with the control signals provided by central plant controller102. For example, BAS 108 may operate equipment 60 using closed loopcontrol to achieve the setpoints specified by central plant controller102. In various embodiments, BAS 108 may be combined with central plantcontroller 102 or may be part of a separate building management system.According to an exemplary embodiment, BAS 108 is a METASYS® brandbuilding management system, as sold by Johnson Controls, Inc.

Central plant controller 102 may monitor the status of the controlledbuilding using information received from BAS 108. Central plantcontroller 102 may be configured to predict the thermal energy loads(e.g., heating loads, cooling loads, etc.) of the building for pluralityof time steps in a prediction window (e.g., using weather forecasts froma weather service). Central plant controller 102 may generate on/offdecisions and/or setpoints for equipment 60 to minimize the cost ofenergy consumed by subplants 12-22 to serve the predicted heating and/orcooling loads for the duration of the prediction window. Central plantcontroller 102 may be configured to carry out process 1100 (FIG. 11 )and other processes described herein. According to an exemplaryembodiment, central plant controller 102 is integrated within a singlecomputer (e.g., one server, one housing, etc.). In various otherexemplary embodiments, central plant controller 102 can be distributedacross multiple servers or computers (e.g., that can exist indistributed locations). In another exemplary embodiment, central plantcontroller 102 may be integrated with a smart building manager thatmanages multiple building systems and/or combined with BAS 108.

Central plant controller 102 is shown to include a communicationsinterface 104 and a processing circuit 106. Communications interface 104may include wired or wireless interfaces (e.g., jacks, antennas,transmitters, receivers, transceivers, wire terminals, etc.) forconducting data communications with various systems, devices, ornetworks. For example, communications interface 104 may include anEthernet card and port for sending and receiving data via anEthernet-based communications network and/or a WiFi transceiver forcommunicating via a wireless communications network. Communicationsinterface 104 may be configured to communicate via local area networksor wide area networks (e.g., the Internet, a building WAN, etc.) and mayuse a variety of communications protocols (e.g., BACnet, IP, LON, etc.).

Communications interface 104 may be a network interface configured tofacilitate electronic data communications between central plantcontroller 102 and various external systems or devices (e.g., BAS 108,subplants 12-22, etc.). For example, central plant controller 102 mayreceive information from BAS 108 indicating one or more measured statesof the controlled building (e.g., temperature, humidity, electric loads,etc.) and one or more states of subplants 12-22 (e.g., equipment status,power consumption, equipment availability, etc.). Communicationsinterface 104 may receive inputs from BAS 108 and/or subplants 12-22 andmay provide operating parameters (e.g., on/off decisions, setpoints,etc.) to subplants 12-22 via BAS 108. The operating parameters may causesubplants 12-22 to activate, deactivate, or adjust a setpoint forvarious devices of equipment 60.

Still referring to FIG. 2 , processing circuit 106 is shown to include aprocessor 110 and memory 112. Processor 110 may be a general purpose orspecific purpose processor, an application specific integrated circuit(ASIC), one or more field programmable gate arrays (FPGAs), a group ofprocessing components, or other suitable processing components.Processor 110 may be configured to execute computer code or instructionsstored in memory 112 or received from other computer readable media(e.g., CDROM, network storage, a remote server, etc.).

Memory 112 may include one or more devices (e.g., memory units, memorydevices, storage devices, etc.) for storing data and/or computer codefor completing and/or facilitating the various processes described inthe present disclosure. Memory 112 may include random access memory(RAM), read-only memory (ROM), hard drive storage, temporary storage,non-volatile memory, flash memory, optical memory, or any other suitablememory for storing software objects and/or computer instructions. Memory112 may include database components, object code components, scriptcomponents, or any other type of information structure for supportingthe various activities and information structures described in thepresent disclosure. Memory 112 may be communicably connected toprocessor 110 via processing circuit 106 and may include computer codefor executing (e.g., by processor 106) one or more processes describedherein.

Still referring to FIG. 2 , memory 112 is shown to include a buildingstatus monitor 134. Central plant controller 102 may receive dataregarding the overall building or building space to be heated or cooledwith central plant 10 via building status monitor 134. In an exemplaryembodiment, building status monitor 134 may include a graphical userinterface component configured to provide graphical user interfaces to auser for selecting building requirements (e.g., overall temperatureparameters, selecting schedules for the building, selecting differenttemperature levels for different building zones, etc.).

Central plant controller 102 may determine on/off configurations andoperating setpoints to satisfy the building requirements received frombuilding status monitor 134. In some embodiments, building statusmonitor 134 receives, collects, stores, and/or transmits cooling loadrequirements, building temperature setpoints, occupancy data, weatherdata, energy data, schedule data, and other building parameters. In someembodiments, building status monitor 134 stores data regarding energycosts, such as pricing information available from utilities 126 (energycharge, demand charge, etc.).

Still referring to FIG. 2 , memory 112 is shown to include a load/rateprediction module 122. Load/rate prediction module 122 may be configuredto predict the thermal energy loads (

) of the building or campus for each time step k (e.g., k=1 . . . n) ofan optimization period. Load/rate prediction module 122 is shownreceiving weather forecasts from a weather service 124. In someembodiments, load/rate prediction module 122 predicts the thermal energyloads

as a function of the weather forecasts. In some embodiments, load/rateprediction module 122 uses feedback from BAS 108 to predict loads

. Feedback from BAS 108 may include various types of sensory inputs(e.g., temperature, flow, humidity, enthalpy, etc.) or other datarelating to the controlled building (e.g., inputs from a HVAC system, alighting control system, a security system, a water system, etc.).

In some embodiments, load/rate prediction module 122 receives a measuredelectric load and/or previous measured load data from BAS 108 (e.g., viabuilding status monitor 134). Load/rate prediction module 122 maypredict loads

as a function of a given weather forecast ({circumflex over (ϕ)}_(w)), aday type (day), the time of day (t), and previous measured load data(Y_(k-1)). Such a relationship is expressed in the following equation:

=f({circumflex over (ϕ)}_(w),day,t|Y _(k-1))

In some embodiments, load/rate prediction module 122 uses adeterministic plus stochastic model trained from historical load data topredict loads

. Load/rate prediction module 122 may use any of a variety of predictionmethods to predict loads

(e.g., linear regression for the deterministic portion and an AR modelfor the stochastic portion). Load/rate prediction module 122 may predictone or more different types of loads for the building or campus. Forexample, load/rate prediction module 122 may predict a hot water load

_(Hot,k) and a cold water load

_(Cold,k) for each time step k within the prediction window.

Load/rate prediction module 122 is shown receiving utility rates fromutilities 126. Utility rates may indicate a cost or price per unit of aresource (e.g., electricity, natural gas, water, etc.) provided byutilities 126 at each time step k in the prediction window. In someembodiments, the utility rates are time-variable rates. For example, theprice of electricity may be higher at certain times of day or days ofthe week (e.g., during high demand periods) and lower at other times ofday or days of the week (e.g., during low demand periods). The utilityrates may define various time periods and a cost per unit of a resourceduring each time period. Utility rates may be actual rates received fromutilities 126 or predicted utility rates estimated by load/rateprediction module 122.

In some embodiments, the utility rates include demand charges for one ormore resources provided by utilities 126. A demand charge may define aseparate cost imposed by utilities 126 based on the maximum usage of aparticular resource (e.g., maximum energy consumption) during a demandcharge period. The utility rates may define various demand chargeperiods and one or more demand charges associated with each demandcharge period. In some instances, demand charge periods may overlappartially or completely with each other and/or with the predictionwindow. Advantageously, optimization module 128 may be configured toaccount for demand charges in the high level optimization processperformed by high level optimization module 130. Utilities 126 may bedefined by time-variable (e.g., hourly) prices, a maximum service level(e.g., a maximum rate of consumption allowed by the physicalinfrastructure or by contract) and, in the case of electricity, a demandcharge or a charge for the peak rate of consumption within a certainperiod.

Load/rate prediction module 122 may store the predicted loads

and the utility rates in memory 112 and/or provide the predicted loads

and the utility rates to optimization module 128. Optimization module128 may use the predicted loads

and the utility rates to determine an optimal load distribution forsubplants 12-22 and to generate on/off decisions and setpoints forequipment 60.

Still referring to FIG. 2 , memory 112 is shown to include anoptimization module 128.

Optimization module 128 may perform a cascaded optimization process tooptimize the performance of central plant 10. For example, optimizationmodule 128 is shown to include a high level optimization module 130 anda low level optimization module 132. High level optimization module 130may control an outer (e.g., subplant level) loop of the cascadedoptimization. High level optimization module 130 may determine anoptimal distribution of thermal energy loads across subplants 12-22 foreach time step in the prediction window in order to optimize (e.g.,minimize) the cost of energy consumed by subplants 12-22. Low leveloptimization module 132 may control an inner (e.g., equipment level)loop of the cascaded optimization. Low level optimization module 132 maydetermine how to best run each subplant at the load setpoint determinedby high level optimization module 130. For example, low leveloptimization module 132 may determine on/off states and/or operatingsetpoints for various devices of equipment 60 in order to optimize(e.g., minimize) the energy consumption of each subplant while meetingthe thermal energy load setpoint for the subplant. The cascadedoptimization process is described in greater detail with reference toFIG. 3 .

Still referring to FIG. 2 , memory 112 is shown to include a subplantcontrol module 138. Subplant control module 138 may store historicaldata regarding past operating statuses, past operating setpoints, andinstructions for calculating and/or implementing control parameters forsubplants 12-22. Subplant control module 138 may also receive, store,and/or transmit data regarding the conditions of individual devices ofequipment 60, such as operating efficiency, equipment degradation, adate since last service, a lifespan parameter, a condition grade, orother device-specific data. Subplant control module 138 may receive datafrom subplants 12-22 and/or BAS 108 via communications interface 104.Subplant control module 138 may also receive and store on/off statusesand operating setpoints from low level optimization module 132.

Data and processing results from optimization module 128, subplantcontrol module 138, or other modules of central plant controller 102 maybe accessed by (or pushed to) monitoring and reporting applications 136.Monitoring and reporting applications 136 may be configured to generatereal time “system health” dashboards that can be viewed and navigated bya user (e.g., a central plant engineer). For example, monitoring andreporting applications 136 may include a web-based monitoringapplication with several graphical user interface (GUI) elements (e.g.,widgets, dashboard controls, windows, etc.) for displaying keyperformance indicators (KPI) or other information to users of a GUI. Inaddition, the GUI elements may summarize relative energy use andintensity across central plants in different buildings (real ormodeled), different campuses, or the like. Other GUI elements or reportsmay be generated and shown based on available data that allow users toassess performance across one or more central plants from one screen.The user interface or report (or underlying data engine) may beconfigured to aggregate and categorize operating conditions by building,building type, equipment type, and the like. The GUI elements mayinclude charts or histograms that allow the user to visually analyze theoperating parameters and power consumption for the devices of thecentral plant.

Still referring to FIG. 2 , central plant controller 102 may include oneor more GUI servers, web services 114, or GUI engines 116 to supportmonitoring and reporting applications 136. In various embodiments,applications 136, web services 114, and GUI engine 116 may be providedas separate components outside of central plant controller 102 (e.g., aspart of a smart building manager). Central plant controller 102 may beconfigured to maintain detailed historical databases (e.g., relationaldatabases, XML databases, etc.) of relevant data and includes computercode modules that continuously, frequently, or infrequently query,aggregate, transform, search, or otherwise process the data maintainedin the detailed databases. Central plant controller 102 may beconfigured to provide the results of any such processing to otherdatabases, tables, XML files, or other data structures for furtherquerying, calculation, or access by, for example, external monitoringand reporting applications.

Central plant controller 102 is shown to include configuration tools118. Configuration tools 118 can allow a user to define (e.g., viagraphical user interfaces, via prompt-driven “wizards,” etc.) howcentral plant controller 102 should react to changing conditions in thecentral plant subsystems. In an exemplary embodiment, configurationtools 118 allow a user to build and store condition-response scenariosthat can cross multiple central plant devices, multiple buildingsystems, and multiple enterprise control applications (e.g., work ordermanagement system applications, entity resource planning applications,etc.). For example, configuration tools 118 can provide the user withthe ability to combine data (e.g., from subsystems, from eventhistories) using a variety of conditional logic. In varying exemplaryembodiments, the conditional logic can range from simple logicaloperators between conditions (e.g., AND, OR, XOR, etc.) to pseudo-codeconstructs or complex programming language functions (allowing for morecomplex interactions, conditional statements, loops, etc.).Configuration tools 118 can present user interfaces for building suchconditional logic. The user interfaces may allow users to definepolicies and responses graphically. In some embodiments, the userinterfaces may allow a user to select a pre-stored or pre-constructedpolicy and adapt it or enable it for use with their system.

Referring now to FIG. 3 , a block diagram illustrating a portion ofcentral plant system 100 in greater detail is shown, according to anexemplary embodiment. FIG. 3 illustrates the cascaded optimizationprocess performed by optimization module 128 to optimize the performanceof central plant 10. In the cascaded optimization process, high leveloptimization module 130 performs a subplant level optimization thatdetermines an optimal distribution of thermal energy loads acrosssubplants 12-22 for each time step in the prediction window in order tominimize the cost of energy consumed by subplants 12-22. Low leveloptimization module 132 performs an equipment level optimization thatdetermines how to best run each subplant at the subplant load setpointdetermined by high level optimization module 130. For example, low leveloptimization module 132 may determine on/off states and/or operatingsetpoints for various devices of equipment 60 in order to optimize theenergy consumption of each subplant while meeting the thermal energyload setpoint for the subplant.

One advantage of the cascaded optimization process performed byoptimization module 128 is the optimal use of computational time. Forexample, the subplant level optimization performed by high leveloptimization module 130 may use a relatively long time horizon due tothe operation of the thermal energy storage. However, the equipmentlevel optimization performed by low level optimization module 132 mayuse a much shorter time horizon or no time horizon at all since the lowlevel system dynamics are relatively fast (compared to the dynamics ofthe thermal energy storage) and the low level control of equipment 60may be handled by BAS 108. Such an optimal use of computational timemakes it possible for optimization module 128 to perform the centralplant optimization in a short amount of time, allowing for real-timepredictive control. For example, the short computational time enablesoptimization module 128 to be implemented in a real-time planning toolwith interactive feedback.

Another advantage of the cascaded optimization performed by optimizationmodule 128 is that the central plant optimization problem can be splitinto two cascaded subproblems. The cascaded configuration provides alayer of abstraction that allows high level optimization module 130 todistribute the thermal energy loads across subplants 12-22 withoutrequiring high level optimization module 130 to know or use any detailsregarding the particular equipment configuration within each subplant.The interconnections between equipment 60 within each subplant may behidden from high level optimization module 130 and handled by low leveloptimization module 132. For purposes of the subplant level optimizationperformed by high level optimization module 130, each subplant may becompletely defined by one or more subplant curves 140.

Still referring to FIG. 3 , low level optimization module 132 maygenerate and provide subplant curves 140 to high level optimizationmodule 130. Subplant curves 140 may indicate the rate of utility use byeach of subplants 12-22 (e.g., electricity use measured in kW, water usemeasured in L/s, etc.) as a function of the subplant load. Exemplarysubplant curves are shown and described in greater detail with referenceto FIGS. 5A-8 . In some embodiments, low level optimization module 132generates subplant curves 140 based on equipment models 120 (e.g., bycombining equipment models 120 for individual devices into an aggregatecurve for the subplant). Low level optimization module 132 may generatesubplant curves 140 by running the low level optimization process forseveral different loads and weather conditions to generate multiple datapoints. Low level optimization module 132 may fit a curve to the datapoints to generate subplant curves 140. In other embodiments, low leveloptimization module 132 provides the data points to high leveloptimization module 132 and high level optimization module 132 generatesthe subplant curves using the data points.

High level optimization module 130 may receive the load and ratepredictions from load/rate prediction module 122 and the subplant curves140 from low level optimization module 132. The load predictions may bebased on weather forecasts from weather service 124 and/or informationfrom building automation system 108 (e.g., a current electric load ofthe building, measurements from the building, a history of previousloads, a setpoint trajectory, etc.). The utility rate predictions may bebased on utility rates received from utilities 126 and/or utility pricesfrom another data source. High level optimization module 130 maydetermine the optimal load distribution for subplants 12-22 (e.g., asubplant load for each subplant) for each time step the predictionwindow and provide the subplant loads as setpoints to low leveloptimization module 132. In some embodiments, high level optimizationmodule 130 determines the subplant loads by minimizing the totaloperating cost of central plant 10 over the prediction window. In otherwords, given a predicted load and utility rate information fromload/rate prediction module 122, high level optimization module 130 maydistribute the predicted load across subplants 12-22 over theoptimization period to minimize operating cost.

In some instances, the optimal load distribution may include using TESsubplants 20 and/or 22 to store thermal energy during a first time stepfor use during a later time step. Thermal energy storage mayadvantageously allow thermal energy to be produced and stored during afirst time period when energy prices are relatively low and subsequentlyretrieved and used during a second time period when energy proves arerelatively high. The high level optimization may be different from thelow level optimization in that the high level optimization has a longertime constant due to the thermal energy storage provided by TESsubplants 20-22. The high level optimization may be described by thefollowing equation:

$\theta_{HL}^{*} = {\arg\min\limits_{\theta_{HL}}{J_{HL}\left( \theta_{HL} \right)}}$where θ_(HL)* contains the optimal high level decisions (e.g., theoptimal load for each of subplants 12-22) for the entire optimizationperiod and J_(HL) is the high level cost function.

To find the optimal high level decisions θ_(HL)*, high leveloptimization module 132 may minimize the high level cost functionJ_(HL). The high level cost function J_(HL) may be the sum of theeconomic costs of each utility consumed by each of subplants 12-22 forthe duration of the optimization period. In some embodiments, the highlevel cost function J_(HL) may be described using the followingequation:

${J_{HL}\left( \theta_{HL} \right)} = {\sum\limits_{k = 1}^{n_{h}}{\sum\limits_{i = 1}^{n_{s}}\left\lbrack {\sum\limits_{j = 1}^{n_{u}}{{t_{s} \cdot c_{jk}}{u_{jik}\left( \theta_{HL} \right)}}} \right\rbrack}}$where n_(h) is the number of time steps k in the optimization period,n_(s) is the number of subplants, t_(s) is the duration of a time step,c_(jk) is the economic cost of utility j at a time step k of theoptimization period, and u_(jik) is the rate of use of utility j bysubplant i at time step k.

In some embodiments, the cost function J_(HL) includes an additionaldemand charge term such as:

$w_{d}c_{demand}{\max\limits_{n_{h}}\left( {{u_{elec}\left( \theta_{HL} \right)},u_{\max,{ele}}} \right)}$where w_(d) is a weighting term, c_(demand) is the demand cost, and themax( ) term selects the peak electricity use during the applicabledemand charge period. Accordingly, the high level cost function J_(HL)may be described by the equation:

${J_{HL}\left( \theta_{HL} \right)} = {{\sum\limits_{k = 1}^{n_{h}}{\sum\limits_{i = 1}^{n_{s}}\left\lbrack {\sum\limits_{j = 1}^{n_{u}}{{t_{s} \cdot c_{jk}}{u_{jik}\left( \theta_{HL} \right)}}} \right\rbrack}} + {w_{d}c_{demand}{\max\limits_{n_{h}}\left( {{u_{elec}\left( \theta_{HL} \right)},u_{\max,{ele}}} \right)}}}$

The decision vector θ_(HL) may be subject to several constraints. Forexample, the constraints may require that the subplants not operate atmore than their total capacity, that the thermal storage not charge ordischarge too quickly or under/over flow for the tank, and that thethermal energy loads for the building or campus are met. Theserestrictions lead to both equality and inequality constraints on thehigh level optimization problem, as described in greater detail withreference to FIG. 4 .

Still referring to FIG. 3 , low level optimization module 132 may usethe subplant loads determined by high level optimization module 130 todetermine optimal low level decisions θ_(LL)* (e.g. binary on/offdecisions, flow setpoints, temperature setpoints, etc.) for equipment60. The low level optimization process may be performed for each ofsubplants 12-22. Low level optimization module 132 may be responsiblefor determining which devices of each subplant to use and/or theoperating setpoints for such devices that will achieve the subplant loadsetpoint while minimizing energy consumption. The low level optimizationmay be described using the following equation:

$\theta_{LL}^{*} = {\arg{\min\limits_{\theta_{LL}}{J_{LL}\left( \theta_{LL} \right)}}}$where θ_(LL)* contains the optimal low level decisions and J_(LL) is thelow level cost function.

To find the optimal low level decisions θ_(LL)*, low level optimizationmodule 132 may minimize the low level cost function J_(LL). The lowlevel cost function J_(LL) may represent the total energy consumptionfor all of equipment 60 in the applicable subplant. The low level costfunction J_(LL) may be described using the following equation:

${J_{LL}\left( \theta_{LL} \right)} = {\sum\limits_{j = 1}^{N}{t_{s} \cdot b_{j} \cdot {u_{j}\left( \theta_{LL} \right)}}}$where N is the number of devices of equipment 60 in the subplant, t_(s)is the duration of a time step, b_(j) is a binary on/off decision (e.g.,0=off, 1=on), and u_(j) is the energy used by device j as a function ofthe setpoint θ_(LL). Each device may have continuous variables which canbe changed to determine the lowest possible energy consumption for theoverall input conditions.

Low level optimization module 132 may minimize the low level costfunction J_(LL) subject to inequality constraints based on thecapacities of equipment 60 and equality constraints based on energy andmass balances. In some embodiments, the optimal low level decisionsθ_(LL)* are constrained by switching constraints defining a shorthorizon for maintaining a device in an on or off state after a binaryon/off switch. The switching constraints may prevent devices from beingrapidly cycled on and off. In some embodiments, low level optimizationmodule 132 performs the equipment level optimization without consideringsystem dynamics. The optimization process may be slow enough to safelyassume that the equipment control has reached its steady-state. Thus,low level optimization module 132 may determine the optimal low leveldecisions θ_(LL)* at an instance of time rather than over a longhorizon.

Low level optimization module 132 may determine optimum operatingstatuses (e.g., on or off) for a plurality of devices of equipment 60.According to an exemplary embodiment, the on/off combinations may bedetermined using binary optimization and quadratic compensation. Binaryoptimization may minimize a cost function representing the powerconsumption of devices in the applicable subplant. In some embodiments,non-exhaustive (i.e., not all potential combinations of devices areconsidered) binary optimization is used. Quadratic compensation may beused in considering devices whose power consumption is quadratic (andnot linear). Low level optimization module 132 may also determineoptimum operating setpoints for equipment using nonlinear optimization.Nonlinear optimization may identify operating setpoints that furtherminimize the low level cost function J_(LL). Low level optimizationmodule 132 may provide the on/off decisions and setpoints to buildingautomation system 108 for use in controlling the central plant equipment60.

In some embodiments, the low level optimization performed by low leveloptimization module 132 is the same or similar to the low leveloptimization process described in U.S. patent application Ser. No.14/634,615 titled “Low Level Central Plant Optimization” and filed onFeb. 27, 2015. The entire disclosure of U.S. patent application Ser. No.14/634,615 is incorporated by reference herein.

Referring now to FIG. 4 , a block diagram illustrating high leveloptimization module 130 in greater detail is shown, according to anexemplary embodiment. High level optimization module 130 may receiveload and rate predictions from load/rate prediction module 122 andsubplant curves from low level optimization module 132. High leveloptimization module 130 may determine optimal subplant loads for each ofsubplants 12-22 as a function of the load and rate predictions and thesubplant curves. In some embodiments, the optimal subplant loadsminimize the economic cost of operating subplants 12-22 to satisfy thepredicted loads for the building or campus. High level optimizationmodule 130 may output the optimal subplant loads to low leveloptimization module 132.

High level optimization module 130 is shown to include an optimizationframework module 142. Optimization framework module 142 may beconfigured to select and/or establish an optimization framework for usein calculating the optimal subplant loads. In some embodiments,optimization framework module 142 uses linear programming as theoptimization framework. A linear programming problem has the followingform:

${{\underset{x}{\arg\min}c^{T}x};{{{subject}{to}{Ax}} \leq b}},{{Hx} = g}$where c is a cost vector, x is a decision matrix, A and b are a matrixand vector (respectively) which describe inequality constraints on theoptimization problem, and H and g are a matrix and vector (respectively)which describe equality constraints on the optimization problem.

The following paragraphs describe an exemplary linear optimizationframework that may be used by high level optimization module 130 tocalculate the optimal subplant loads. Advantageously, the linearprogramming framework described herein allows high level optimizationmodule 130 to determine the subplant load distribution for a longoptimization period in a very short timeframe complete with load changepenalties, demand charges, and subplant performance curves. However, thelinear optimization framework is merely one example of an optimizationframework that can be used by high level optimization module 130 andshould not be regarded as limiting. It should be understood that inother embodiments, high level optimization module 130 may use any of avariety of other optimization frameworks and/or optimization techniques(e.g., quadratic programming, linear-fractional programming, nonlinearprogramming, combinatorial algorithms, etc.) to calculate the optimalsubplant loads.

Still referring to FIG. 4 , high level optimization module 130 is shownto include a linear program module 144. Linear program module 144 may beconfigured to formulate and solve a linear optimization problem tocalculate the optimal subplant loads. For example, linear program module144 may determine and set values for the cost vector c, the A matrix andthe b vector which describe the inequality constraints, and the H matrixand the g vector which describe the equality constraints. Linear programmodule 144 may determine an optimal decision matrix x* that minimizesthe cost function c^(T)x. The optimal decision matrix x* may correspondto the optimal decisions θ_(HL)* (for each time step k within anoptimization period) that minimize the high level cost function J_(HL),as described with reference to FIG. 3 .

For a central plant 10 that includes chillers, heat recovery chillers,hot water generators, and thermal energy storage, the plant assetsacross which the loads are to be distributed may include a chillersubplant 16, a heat recovery chiller subplant 14, a heater subplant 12,a hot thermal energy storage subplant 20, and a cold thermal energystorage subplant 22. The loads across each of subplants 12-22 may be thedecision variables in the decision matrix x that the high leveloptimization determines for each time step k within the optimizationperiod. For example, linear program module 144 may formulate thedecision matrix x as:x=[{dot over (Q)} _(Chiller,1 . . . n) , {dot over (Q)}_(hrChiller,1 . . . n) , {dot over (Q)} _(Heater,1 . . . n) , {dot over(Q)} _(HotStorage,1 . . . n) , {dot over (Q)}_(ColdStorage,1 . . . n)]^(T)where {dot over (Q)}_(Chiller,1 . . . n), {dot over(Q)}_(hrChiller,1 . . . n), {dot over (Q)}_(Heater,1 . . . n), {dot over(Q)}_(HotStorage,1, . . . n), and {dot over (Q)}_(ColdStorage,1 . . . n)are n-dimensional vectors representing the thermal energy load assignedto chiller subplant 16, heat recovery chiller subplant 14, heatersubplant 12, hot TES subplant 20, and cold TES subplant 22,respectively, for each of the n time steps within the optimizationperiod.

Linear program module 144 may formulate the linear program for thesimple case where only energy cost and equipment constraints areconsidered. The simplified linear program may then be modified byinequality constraints module 146, equality constraints module 148,unmet loads module 150, ground loop module 152, heat exchanger module154, demand charge module 156, load change penalty module 158, tankforced full module 160, and/or subplant curves module 170 to provideadditional enhancements, described in greater detail below.

In some embodiments, linear program module 144 formulates the simplifiedlinear program using the assumption that each subplant has a specificcost per unit load. For example, linear program module 144 may assumethat each subplant has a constant coefficient of performance (COP) orefficiency for any given time step k. The COP can change over time andmay have a different value for different time steps; however, in thesimplest case, the COP for each of subplant is not a function of theloading. With this assumption, linear program module 144 may formulatethe cost function c as:

$c = {t_{s} \cdot \left\lbrack {\left\lbrack {\sum\limits_{j = 1}^{n_{u}}{c_{j}u_{j,{Chiller}}}} \right\rbrack_{1\ldots h},\left\lbrack {\sum\limits_{j = 1}^{n_{u}}{c_{j}u_{j,{hrChiller}}}} \right\rbrack_{1\ldots h},\left\lbrack {\sum\limits_{j = 1}^{n_{u}}{c_{j}u_{j,{Heater}}}} \right\rbrack_{1\ldots h},0_{h},0_{h}} \right\rbrack^{T}}$where t_(s) is the duration of a time step, n_(u) is the total number ofresources (e.g., electricity, natural gas, water, etc.) consumed by thesubplants, c_(j) is the cost per unit of the jth resource, andu_(j,Chiller), u_(j,hrChiller), and u_(j,Heater) are the usage rates ofthe jth resource by chiller subplant 16, heat recovery chiller subplant14, and heater subplant 12, respectively, for each of the h time stepswithin the optimization period. The first three elements of the form[Σ_(j=1) ^(n) ^(u) c_(j)u_(j)]_(1 . . . h) represent vectors of h sums(i.e., summing over all resource use), one for each time step within theoptimization period. The last two elements of the form 0_(h) are zero toindicate that charging or discharging the thermal energy storage tankshas no cost (pumping power is neglected).

In some embodiments, linear program module 144 uses the load and ratepredictions to formulate the linear program. For example, linear programmodule 144 may use the load predictions to determine values foru_(j,Chiller), u_(j,hrChiller), and u_(j,Heater) and may use the ratepredictions to determine values for c_(j) for each of the n_(u)resources. In some embodiments, linear program module 144 uses thesubplant curves to define c_(j) as a function of the resource usage.Linear program module 144 may use inputs from inequality constraintsmodule 146, equality constraints module 148, unmet loads module 150,ground loop module 152, heat exchanger module 154, demand charge module156, load change penalty module 158, tank forced full module 160, and/orsubplant curves module 170 to determine and set values for the variousmatrices and vectors in the linear program. Modules 146-170 may modifythe cost vector c, the A matrix, the b vector, the H matrix, and/or theg vector to provide additional enhancements and/or functionality to thelinear program. The inputs provided by modules 146-170 are described ingreater detail below.

Linear program module 144 may use any of a variety of linearoptimization techniques to solve the linear optimization problem. Forexample, linear program module 144 may use basis exchange algorithms(e.g., simplex, crisscross, etc.), interior point algorithms (e.g.,ellipsoid, projective, path-following, etc.), covering and packingalgorithms, integer programming algorithms (e.g., cutting-plant, branchand bound, branch and cut, branch and price, etc.), or any other type oflinear optimization algorithm or technique to solve the linear programsubject to the optimization constraints. For embodiments in whichnonlinear optimization is used, linear program module 144 may use any ofa variety of nonlinear optimization techniques to solve the nonlinearoptimization problem.

Still referring to FIG. 4 , high level optimization module 130 is shownto include an inequality constraints module 146. Inequality constraintsmodule 146 may formulate or define one or more inequality constraints onthe optimization problem solved by linear program module 144. In someinstances, inequality constraints module 146 defines inequalityconstraints on the decision variables {dot over (Q)}_(Chiller,k), {dotover (Q)}_(hrChiller,k), and {dot over (Q)}_(Heater,k) corresponding tothe loads on chiller subplant 16, heat recovery chiller subplant 14, andheater subplant 12, respectively, for each time step k withinoptimization period. For example, each of subplants 12-16 may have twocapacity constraints given by the following equations:{dot over (Q)} _(i,k) ≤{dot over (Q)} _(i,max) ∀k∈horizon{dot over (Q)} _(i,k)≥0∀k∈horizonwhere {dot over (Q)}_(i,k) is the load on the ith subplant during timestep k and {dot over (Q)}_(i,max) is the maximum capacity of the ithsubplant. The first capacity constraint requires the load {dot over(Q)}_(i,k) on each of subplants 12-16 to be less than or equal to themaximum capacity i,max of the subplant for each time step k within theoptimization period. The second capacity constraint requires the load{dot over (Q)}_(i,k) on each of subplants 12-16 to be greater than orequal to zero for each time step k within the optimization period.

The inequality constraints for chiller subplant 16 can be placed in theform Ax≤b by defining the A matrix and the b vector as follows:

${A = \begin{bmatrix}\left\lbrack I_{h} \right\rbrack & \left\lbrack 0_{h} \right\rbrack & \left\lbrack 0_{h} \right\rbrack & \left\lbrack 0_{h} \right\rbrack & \left\lbrack 0_{h} \right\rbrack \\\left\lbrack {- I_{h}} \right\rbrack & \left\lbrack 0_{h} \right\rbrack & \left\lbrack 0_{h} \right\rbrack & \left\lbrack 0_{h} \right\rbrack & \left\lbrack 0_{h} \right\rbrack\end{bmatrix}},{b = \begin{bmatrix}{{\overset{.}{Q}}_{{Chiller},\max}\left\lbrack I_{h} \right\rbrack} \\\left\lbrack 0_{h} \right\rbrack\end{bmatrix}}$where [I_(h)] represents either an h by h identity matrix or an h by 1ones vector, [0_(h)] represents either an h by h zero matrix or an h by1 zero vector, and {dot over (Q)}_(Chiller,max) is the maximum capacityof chiller subplant 16. Similar inequality constraints for heat recoverychiller subplant 14 and heater subplant 12 can be placed in the formAx≤b by defining the A matrices and the b vectors as follows:

${{A = \begin{bmatrix}\left\lbrack 0_{h} \right\rbrack & \left\lbrack I_{h} \right\rbrack & \left\lbrack 0_{h} \right\rbrack & \left\lbrack 0_{h} \right\rbrack & \left\lbrack 0_{h} \right\rbrack \\\left\lbrack 0_{h} \right\rbrack & \left\lbrack {- I_{h}} \right\rbrack & \left\lbrack 0_{h} \right\rbrack & \left\lbrack 0_{h} \right\rbrack & \left\lbrack 0_{h} \right\rbrack\end{bmatrix}},{b = \begin{bmatrix}{{\overset{.}{Q}}_{{hrChiller},\max}\left\lbrack I_{h} \right\rbrack} \\\left\lbrack 0_{h} \right\rbrack\end{bmatrix}}}{{A = \begin{bmatrix}\left\lbrack 0_{h} \right\rbrack & \left\lbrack 0_{h} \right\rbrack & \left\lbrack I_{h} \right\rbrack & \left\lbrack 0_{h} \right\rbrack & \left\lbrack 0_{h} \right\rbrack \\\left\lbrack 0_{h} \right\rbrack & \left\lbrack 0_{h} \right\rbrack & \left\lbrack {- I_{h}} \right\rbrack & \left\lbrack 0_{h} \right\rbrack & \left\lbrack 0_{h} \right\rbrack\end{bmatrix}},{b = \begin{bmatrix}{{\overset{.}{Q}}_{{Heater},\max}\left\lbrack I_{h} \right\rbrack} \\\left\lbrack 0_{h} \right\rbrack\end{bmatrix}}}$where {dot over (Q)}_(hrChiller,max) is the maximum capacity of heatrecovery chiller subplant 14 and {dot over (Q)}_(Heater,max) is themaximum capacity of heater subplant 12.

Inequality constraints module 146 may formulate or define inequalityconstraints on the decision variables {dot over (Q)}_(HotStorage,k) and{dot over (Q)}_(ColdStorage,k) corresponding to the loads on hot TESsubplant 20 and cold TES subplant 22 for each time step k within theoptimization period. For example, each of subplants 20-22 may have twocapacity constraints given by the following equations:{dot over (Q)} _(i,k) ≤{dot over (Q)} _(discharge,i,max) ∀k∈horizon−{dot over (Q)} _(i,k) ≤{dot over (Q)} _(charge,i,max) ∀k∈horizonwhere {dot over (Q)}_(i,k) is the rate at which ith TES subplant isbeing discharged at time step k, {dot over (Q)}_(discharge,i,max) is themaximum discharge rate of the ith subplant, and {dot over(Q)}_(charge,i,max) is the maximum charge rate of the ith subplant.Positive load values for {dot over (Q)}_(i,k) indicate that the TESsubplant is discharging and negative load values for {dot over(Q)}_(i,k) indicate that the subplant is charging. The first capacityconstraint requires the discharge rate {dot over (Q)}_(i,k) for each ofsubplants 20-22 to be less than or equal to the maximum discharge rate{dot over (Q)}_(discharge,i,max) of the subplant for each time step kwithin the optimization period. The second capacity constraint requiresthe negative discharge rate −{dot over (Q)}_(i,k) (i.e., the chargerate) for each of subplants 20-22 to be less than or equal to themaximum charge rate {dot over (Q)}_(charge,i,max) of the subplant foreach time step k within the optimization period.

The inequality constraints for hot TES subplant 20 can be placed in theform Ax≤b by defining the A matrix and the b vector as follows:

${A = \begin{bmatrix}\left\lbrack 0_{h} \right\rbrack & \left\lbrack 0_{h} \right\rbrack & \left\lbrack 0_{h} \right\rbrack & \left\lbrack I_{h} \right\rbrack & \left\lbrack 0_{h} \right\rbrack \\\left\lbrack 0_{h} \right\rbrack & \left\lbrack 0_{h} \right\rbrack & \left\lbrack 0_{h} \right\rbrack & \left\lbrack {- I_{h}} \right\rbrack & \left\lbrack 0_{h} \right\rbrack\end{bmatrix}},{b = \begin{bmatrix}{{\overset{.}{Q}}_{{HotDischarge},\max}\left\lbrack I_{h} \right\rbrack} \\{{\overset{.}{Q}}_{{HotCharge},\max}\left\lbrack I_{h} \right\rbrack}\end{bmatrix}}$where {dot over (Q)}_(HotDischarge,max) is the maximum discharge ratefor hot TES subplant 20 and {dot over (Q)}_(HotCharge,max) is themaximum charge rate for hot TES subplant 20. Similar inequalityconstraints for cold TES subplant 22 can be placed in the form Ax≤b bydefining the A matrix and the b vector as follows:

${A = \begin{bmatrix}\left\lbrack 0_{h} \right\rbrack & \left\lbrack 0_{h} \right\rbrack & \left\lbrack 0_{h} \right\rbrack & \left\lbrack 0_{h} \right\rbrack & \left\lbrack I_{h} \right\rbrack \\\left\lbrack 0_{h} \right\rbrack & \left\lbrack 0_{h} \right\rbrack & \left\lbrack 0_{h} \right\rbrack & \left\lbrack 0_{h} \right\rbrack & \left\lbrack {- I_{h}} \right\rbrack\end{bmatrix}},{b = \begin{bmatrix}{{\overset{.}{Q}}_{{ColdDischarge},\max}\left\lbrack I_{h} \right\rbrack} \\{{\overset{.}{Q}}_{{ColdCharge},\max}\left\lbrack I_{h} \right\rbrack}\end{bmatrix}}$where {dot over (Q)}_(ColdDischarge,max) is the maximum discharge ratefor cold TES subplant 22 and {dot over (Q)}_(ColdCharge,max) is themaximum charge rate for cold TES subplant 22.

Inequality constraints module 146 may implement an electrical demandconstraint for the total electrical usage of all the subplants and thebuilding/campus P_(elec,campus). Inequality constraints module 146 mayrequire that the total electrical demand be less than or equal to amaximum electrical demand P_(elec,max) by defining the A matrix and theb vector as follows:A=[u _(elec,Chiller)[I _(h)],u _(elec,hrChiller)[I _(h)],u_(elec,Heater)[I _(h)],0_(n),0_(n)],b=P _(elec,max)[I _(h)]−P_(elec,campus,k)where u_(elec,Chiller), u_(elec,hrChiller), and u_(elec,Heater) are theelectrical usage values for chiller subplant 16, heat recovery chillersubplant 14, and heater subplant 12, respectively, P_(elec,campus,k) isthe electrical usage of the building/campus at time k, and P_(elec,max)is the maximum total electrical usage for central plant 10 and thebuilding/campus.

Inequality constraints module 146 may implement tank capacityconstraints for hot TES subplant 20 and cold TES subplant 22. The tankcapacity constraints may require that each TES tank never charge aboveits maximum capacity or discharge below zero. These physicalrequirements lead to a series of constraints to ensure that the initialtank level Q₀ of each TES tank at the beginning of the optimizationperiod plus all of the charging during time steps 1 to k into theoptimization period is less than or equal to the maximum capacityQ_(max) of the TES tank. A similar constraint may be implemented toensure that the initial tank level Q₀ of each TES tank at the beginningof the optimization period minus all of the discharging during timesteps 1 to k into the optimization period is greater than or equal tozero.

The tank capacity constraints for hot TES subplant 20 can be placed inthe form Ax≤b by defining the A matrix and the b vector as follows:

${A = \begin{bmatrix}\left\lbrack 0_{h} \right\rbrack & \left\lbrack 0_{h} \right\rbrack & \left\lbrack 0_{h} \right\rbrack & {t_{s}\left\lbrack \Delta_{h} \right\rbrack} & \left\lbrack 0_{h} \right\rbrack \\\left\lbrack 0_{h} \right\rbrack & \left\lbrack 0_{h} \right\rbrack & \left\lbrack 0_{h} \right\rbrack & {t_{s}\left\lbrack {- \Delta_{h}} \right\rbrack} & \left\lbrack 0_{h} \right\rbrack\end{bmatrix}},{b = \begin{bmatrix}{Q_{0,{Hot}}\left\lbrack I_{h} \right\rbrack} \\{Q_{\max,{Hot}} - {Q_{0,{Hot}}\left\lbrack I_{h} \right\rbrack}}\end{bmatrix}}$where Q_(0,Hot) is the initial charge level of hot TES subplant 20 atthe beginning of the optimization period, Q_(max,Hot) is the maximumcharge level of hot TES subplant 20, Δ_(h) is a lower triangular matrixof ones, and t_(s) is the duration of a time step. Discharging the tankis represented in the top row of the A matrix as positive flow from thetank and charging the tank is represented in the bottom row of the Amatrix as negative flow from the tank. Similar inequality constraintsfor cold TES subplant 22 can be placed in the form Ax≤b by defining theA matrix and the b vector as follows:

${A = \begin{bmatrix}\left\lbrack 0_{h} \right\rbrack & \left\lbrack 0_{h} \right\rbrack & \left\lbrack 0_{h} \right\rbrack & \left\lbrack 0_{h} \right\rbrack & {t_{s}\left\lbrack \Delta_{h} \right\rbrack} \\\left\lbrack 0_{h} \right\rbrack & \left\lbrack 0_{h} \right\rbrack & \left\lbrack 0_{h} \right\rbrack & \left\lbrack 0_{h} \right\rbrack & {t_{s}\left\lbrack {- \Delta_{h}} \right\rbrack}\end{bmatrix}},{b = \begin{bmatrix}{Q_{0,{Cold}}\left\lbrack I_{h} \right\rbrack} \\{Q_{\max,{Cold}} - {Q_{0,{Cold}}\left\lbrack I_{h} \right\rbrack}}\end{bmatrix}}$where Q_(0,Cold) is the initial charge level of cold TES subplant 22 atthe beginning of the optimization period and Q_(max,Cold) is the maximumcharge level of cold TES subplant 22.

Still referring to FIG. 4 , high level optimization module 130 is shownto include an equality constraints module 148. Equality constraintsmodule 148 may formulate or define one or more equality constraints onthe optimization problem solved by linear program module 144. Theequality constraints may ensure that the predicted thermal energy loadsof the building or campus are satisfied for each time step k in theoptimization period. Equality constraints module 148 may formulate anequality constraint for each type of thermal energy load (e.g., hotwater, cold water, etc.) to ensure that the load is satisfied. Theequality constraints may be given by the following equation:

${\sum\limits_{i = 1}^{n_{s}}{\overset{.}{Q}}_{p,i,k}} = {{\hat{\ell}}_{p,k}{\forall{k \in {horizon}}}}$where {dot over (Q)}_(p,i,k) is the thermal energy load of type p (e.g.,hot water, cold water, etc.) on the ith subplant during time step k,n_(s) is the total number of subplants capable of serving thermal energyload p, and

is the predicted thermal energy load of type p that must be satisfied attime step k. The predicted thermal energy loads may be received as loadpredictions from load/rate prediction module 122.

In some embodiments, the predicted thermal energy loads include apredicted hot water thermal energy load

and a predicted cold water thermal energy load

for each time step k. The predicted hot water thermal energy load

may be satisfied by the combination of heat recovery chiller subplant14, heater subplant 12, and hot TES subplant 20. The predicted coldwater thermal energy load

may be satisfied by the combination of heat recovery chiller subplant14, chiller subplant 16, and cold TES subplant 22.

The equality constraints can be placed in the form Hx=g by defining theH matrix and the g vector as follows:

${H = \begin{bmatrix}\left\lbrack I_{h} \right\rbrack & \left\lbrack I_{h} \right\rbrack & \left\lbrack 0_{h} \right\rbrack & \left\lbrack 0_{h} \right\rbrack & \left\lbrack I_{h} \right\rbrack \\\left\lbrack 0_{h} \right\rbrack & \left( {1 + {u_{{elec},{hrChiller}}\left\lbrack I_{h} \right\rbrack}} \right. & \left\lbrack I_{h} \right\rbrack & \left\lbrack I_{h} \right\rbrack & \left\lbrack 0_{h} \right\rbrack\end{bmatrix}},{g = \begin{bmatrix}{\hat{\ell}}_{{Cold},{1\ldots k}} \\{\hat{\ell}}_{{Hot},{1\ldots k}}\end{bmatrix}}$where

and

are k-dimensional vectors of predicted cold water loads and predictedhot water loads, respectively, at each of time steps k, andu_(elec,hrChiller) is the electrical consumption of heat recoverychiller subplant 14. For central plants that serve one or moreadditional types of loads, an additional row may be added to the Hmatrix and the g vector to define the equality constraints for eachadditional load served by the central plant.

For this example problem, assuming an optimization period of 72 one-hoursamples, the linear program has 360 decision variables and 1224constraints. However, linear program module 144 can solve this linearprogram to determine the optimal subplant load values in less than 200milliseconds using the linear programming framework. Advantageously,this allows high level optimization module 130 to determine the subplantload distribution for a long optimization period in a very shorttimeframe.

Still referring to FIG. 4 , high level optimization module 130 is shownto include an unmet loads module 150. In some instances, the centralplant equipment 60 may not have enough capacity or reserve storage tosatisfy the predicted thermal energy loads, regardless of how thethermal energy loads are distributed across subplants 12-22. In otherwords, the high level optimization problem may have no solution thatsatisfies all of the inequality and equality constraints, even if theapplicable subplants are operated at maximum capacity. Unmet loadsmodule 150 may be configured to modify the high level optimizationproblem to account for this possibility and to allow the high leveloptimization to find the solution that results in the minimal amount ofunmet loads.

In some embodiments, unmet loads module 150 modifies the decisionvariable matrix x by introducing a slack variable for each type ofthermal energy load. The slack variables represent an unsatisfied (e.g.,unmet, deferred, etc.) amount of each type of thermal energy load. Forexample, unmet loads module 150 may modify the decision variable matrixx as follows:x=[{dot over (Q)} _(Chiller,1 . . . n) ,{dot over (Q)}_(hrChiller,1 . . . n) ,{dot over (Q)} _(Heater,1 . . . n) ,{dot over(Q)} _(HotStorage,1 . . . n) ,{dot over (Q)} _(ColdStorage,1 . . . n) ,Q_(ColdUnmet,1 . . . n) ,Q _(HotUnmet,1 . . . n)]^(T)where Q_(ColdUnmet,1 . . . n) and Q_(HotUnmet,1 . . . n) aren-dimensional vectors representing a total deferred cold thermal energyload and a total deferred hot thermal energy load, respectively, at eachtime step k within the optimization period. In some embodiments, thedecision variables Q_(ColdUnmet,1 . . . n) and Q_(HotUnmet,1 . . . n)represent total deferred loads that have accumulated up to each timestep k rather than the incremental deferred load at each time step. Thetotal deferred load may be used because any deferred load is likely toincrease the required load during subsequent time steps.

Unmet loads module 150 may modify the equality constraints to accountfor any deferred thermal energy loads. The modified equality constraintsmay require that the predicted thermal energy loads are equal to thetotal loads satisfied by subplants 12-22 plus any unsatisfied thermalenergy loads. The modified equality constraints can be placed in theform Hx=g by defining the H matrix and the g vector as follows:

$H = {{\begin{bmatrix}\left\lbrack I_{h} \right\rbrack & \left\lbrack I_{h} \right\rbrack & \left\lbrack 0_{h} \right\rbrack & \left\lbrack 0_{h} \right\rbrack & \left\lbrack I_{h} \right\rbrack & {\left\lbrack I_{h} \right\rbrack - \left\lbrack D_{- 1} \right\rbrack} & \left\lbrack 0_{h} \right\rbrack \\\left\lbrack 0_{h} \right\rbrack & \left( {1 + {u_{{elec},{hrChiller}}\left\lbrack I_{h} \right\rbrack}} \right. & \left\lbrack I_{h} \right\rbrack & \left\lbrack I_{h} \right\rbrack & \left\lbrack 0_{h} \right\rbrack & \left\lbrack 0_{h} \right\rbrack & {\left\lbrack I_{h} \right\rbrack - \left\lbrack D_{- 1} \right\rbrack}\end{bmatrix},{g = \begin{bmatrix}{\hat{\ell}}_{{Cold},{1\ldots k}} \\{\hat{\ell}}_{{Hot},{1\ldots k}}\end{bmatrix}}}}$where [D⁻¹] is a lower diagonal matrix of ones.

Unmet loads module 150 may modify the cost vector c to associate costvalues with any unmet loads. In some embodiments, unmet loads module 150assigns unmet loads a relatively higher cost compared to the costsassociated with other types of loads in the decision variable matrix x.Assigning a large cost to unmet loads ensures that the optimal solutionto the high level optimization problem uses unmet loads only as a lastresort (i.e., when the optimization has no solution without using unmetloads). Accordingly, linear program module 144 may avoid using unmetloads if any feasible combination of equipment is capable of satisfyingthe predicted thermal energy loads. In some embodiments, unmet loadsmodule 150 assigns a cost value to unmet loads that allows linearprogram module 144 to use unmet loads in the optimal solution even ifthe central plant is capable of satisfying the predicted thermal energyloads. For example, unmet loads module 150 may assign a cost value thatallows linear program module 144 to use unmet loads if the solutionwithout unmet loads would be prohibitively expensive and/or highlyinefficient.

Still referring to FIG. 4 , high level optimization module 130 is shownto include a subplant curves module 170. In the simplest case describedwith reference to linear program module 144, it was assumed that theresource consumption of each subplant is a linear function of thethermal energy load produced by the subplant. However, this assumptionmay not be true for some subplant equipment, much less for an entiresubplant. Subplant curves module 170 may be configured to modify thehigh level optimization problem to account for subplants that have anonlinear relationship between resource consumption and load production.

Subplant curves module 170 is shown to include a subplant curve updater172, a subplant curves database 174, a subplant curve linearizer 176,and a subplant curves incorporator 178. Subplant curve updater 172 maybe configured to request subplant curves for each of subplants 12-22from low level optimization module 132. Each subplant curve may indicatean amount of resource consumption by a particular subplant (e.g.,electricity use measured in kW, water use measured in L/s, etc.) as afunction of the subplant load. Exemplary subplant curves are shown anddescribed in greater detail with reference to FIGS. 5A-8 .

In some embodiments, low level optimization module 132 generates thesubplant curves by running the low level optimization process forvarious combinations of subplant loads and weather conditions togenerate multiple data points. Low level optimization module 132 may fita curve to the data points to generate the subplant curves and providethe subplant curves to subplant curve updater 172. In other embodiments,low level optimization module 132 provides the data points to subplantcurve updater 172 and subplant curve updater 172 generates the subplantcurves using the data points. Subplant curve updater 172 may store thesubplant curves in subplant curves database 174 for use in the highlevel optimization process.

In some embodiments, the subplant curves are generated by combiningefficiency curves for individual devices of a subplant. A deviceefficiency curve may indicate the amount of resource consumption by thedevice as a function of load. The device efficiency curves may beprovided by a device manufacturer or generated using experimental data.In some embodiments, the device efficiency curves are based on aninitial efficiency curve provided by a device manufacturer and updatedusing experimental data. The device efficiency curves may be stored inequipment models 120. For some devices, the device efficiency curves mayindicate that resource consumption is a U-shaped function of load.Accordingly, when multiple device efficiency curves are combined into asubplant curve for the entire subplant, the resultant subplant curve maybe a wavy curve as shown in FIG. 6 . The waves are caused by a singledevice loading up before it is more efficient to turn on another deviceto satisfy the subplant load.

Subplant curve linearizer 176 may be configured to convert the subplantcurves into convex curves. A convex curve is a curve for which a lineconnecting any two points on the curve is always above or along thecurve (i.e., not below the curve). Convex curves may be advantageous foruse in the high level optimization because they allow for anoptimization process that is less computationally expensive relative toan optimization process that uses non-convex functions. Subplant curvelinearizer 176 may be configured to break the subplant curves intopiecewise linear segments that combine to form a piecewise-definedconvex curve. An unmodified subplant curve 600 and a linearized subplantcurve 700 generated by subplant curve linearizer 176 are shown in FIGS.6 and 7 , respectively. Subplant curve linearizer 176 may store thelinearized subplant curves in subplant curves database 174.

Still referring to FIG. 4 , subplant curves module 170 is shown toinclude a subplant curve incorporator 178. Subplant curve incorporator178 may be configured to modify the high level optimization problem toincorporate the subplant curves into the optimization. In someembodiments, subplant curve incorporator 178 modifies the decisionmatrix x to include one or more decision vectors representing theresource consumption of each subplant. For example, for chiller subplant16, subplant curve incorporator 178 may modify the decision matrix x asfollows:x=[ . . . {dot over (Q)} _(Chiller,1 . . . n) . . . u_(Chiller,elec,1 . . . n) u _(Chiller,water,1 . . . n) . . . ]^(T)where u_(Chiller,elec,1 . . . n) and u_(Chiller,water,1 . . . n) aren-dimensional vectors representing the amount of electrical consumptionand water consumption, respectively, by chiller subplant 16 at each timestep k.

Subplant curve incorporator 178 may add one or more resource consumptionvectors to matrix x for each of subplants 12-22. The decision vectorsadded by subplant curve incorporator 178 for a given subplant mayrepresent an amount of resource consumption for each resource consumedby the subplant (e.g., water, electricity, natural gas, etc.) at eachtime step k within the optimization period. For example, if heatersubplant 12 consumes natural gas, electricity, and water, subplant curveincorporator 178 may add a decision vector u_(Heater,gas,1 . . . n)representing an amount of natural gas consumed by heater subplant 12 ateach time step, a decision vector u_(Heater,elec,1 . . . n) representingan amount of electricity consumed by heater subplant 12 at each timestep, and a decision vector u_(Heater,water,1 . . . n) representing anamount of water consumed by heater subplant at each time step. Subplantcurve incorporator 178 may add resource consumption vectors for othersubplants in a similar manner.

Subplant curve incorporator 178 may modify the cost vector c to accountfor the resource consumption vectors in the decision matrix x. In someembodiments, subplant curve incorporator 178 removes (or sets to zero)any cost directly associated with the subplant loads (e.g., {dot over(Q)}_(Chiller,1 . . . n), {dot over (Q)}_(Heater,1 . . . n), etc.) andadds economic costs associated with the resource consumption required toproduce the subplant loads. For example, for chiller subplant 16,subplant curve incorporator 178 may modify the cost vector c as follows:c=[ . . . 0_(n) . . . c _(elec,1 . . . n) c _(water,1 . . . n) . . .]^(T)where 0_(n) is a n-dimensional zero vector indicating that the directeconomic cost of {dot over (Q)}_(Chiller,1 . . . n) is zero at each timestep, c_(elec,1 . . . n) is a n-dimensional vector indicating the perunit cost of electricity at each time step, and c_(water,1 . . . n) is an-dimensional vector indicating the per unit cost of water at each timestep. The modified cost vector associates an economic cost with theresources consumed to produce the subplant loads rather than thesubplant loads themselves. In some embodiments, the values forc_(elec,1 . . . n) and c_(water,1 . . . n) are utility rates obtainedfrom load/rate prediction module 122.

Subplant curve incorporator 178 may modify the inequality constraints toensure that the proper amount of each resource is consumed to serve thepredicted thermal energy loads. In some embodiments, subplant curveincorporator 178 formulates inequality constraints that force theresource usage for each resource in the epigraph of the correspondinglinearized subplant curve. For example, chiller subplant 16 may have alinearized subplant curve that indicates the electricity use of chillersubplant 16 (i.e., u_(Chiller,elec)) as a function of the cold waterproduction of chiller subplant 16 (i.e., {dot over (Q)}_(Chiller)). Sucha linearized subplant curve 700 is shown in FIG. 7 . The linearizedsubplant curve may include a first line segment connecting point [u₁,Q₁] to point [u₂, Q₂], a second line segment connecting point [u₂, Q₂]to point [u₃, Q₃], and a third line segment connecting point [u₃, Q₃] topoint [u₄, Q₄].

Subplant curve incorporator 178 may formulate an inequality constraintfor each piecewise segment of the subplant curve that constrains thevalue of u_(Chiller,elec) to be greater than or equal to the amount ofelectricity use defined by the line segment for the corresponding valueof {dot over (Q)}_(Chiller). The subplant curve constraints for theelectricity use of chiller subplant 16 can be placed in the form Ax≤b bydefining the A matrix and the b vector as follows:

${A = \begin{bmatrix}\cdots & {\left\lbrack {- \left( {u_{2} - u_{1}} \right)} \right\rbrack I_{n}} & \cdots & {\left\lbrack \left( {Q_{2} - Q_{1}} \right) \right\rbrack I_{n}} & 0_{n} & \cdots \\\cdots & {\left\lbrack {- \left( {u_{3} - u_{2}} \right)} \right\rbrack I_{n}} & \cdots & {\left\lbrack \left( {Q_{3} - Q_{2}} \right) \right\rbrack I_{n}} & 0_{n} & \cdots \\\cdots & {\left\lbrack {- \left( {u_{4} - u_{3}} \right)} \right\rbrack I_{n}} & \cdots & {\left\lbrack \left( {Q_{4} - Q_{3}} \right) \right\rbrack I_{n}} & 0_{n} & \cdots\end{bmatrix}},{b = \begin{bmatrix}{{Q_{1}u_{2}} - {Q_{2}u_{1}}} \\{{Q_{2}u_{3}} - {Q_{3}u_{2}}} \\{{Q_{3}u_{4}} - {Q_{4}u_{3}}}\end{bmatrix}}$

Similar inequality constraints can be formulated for other subplantcurves. For example, subplant curve incorporator 178 may generate a setof inequality constraints for the water consumptionu_(Chiller,water,1 . . . n) of chiller subplant 16 using the pointsdefining the linearized subplant curve for the water consumptionu_(Chiller,water,1 . . . n) of chiller subplant 16 as a function of coldwater production {dot over (Q)}_(Chiller). In some embodiments, thewater consumption of chiller subplant 16 is equal to the cold waterproduction and the linearized subplant curve for water consumptionincludes a single line segment connecting point [u₅, Q₅] to point [u₆,Q₆] (as shown in FIG. 5B). The subplant curve constraints for the coldwater consumption of chiller subplant 16 can be placed in the form Ax≤bby defining the A matrix and the b vector as follows:A=[ . . . [−(u ₆ −u ₅)]I _(n) . . . 0_(n) [(Q ₆ −Q ₅)]I _(n) . . . ],b=[Q ₅ u ₆ −Q ₆ u ₅]Subplant curve incorporator 178 may repeat this process for eachsubplant curve for chiller subplant 16 and for the other subplants ofcentral plant 10 to define a set of inequality constraints for eachsubplant curve.

The inequality constraints generated by subplant curve incorporator 178ensure that high level optimization module 130 keeps the resourceconsumption above all of the line segments of the corresponding subplantcurve. In most situations, there is no reason for high leveloptimization module 130 to choose a resource consumption value that liesabove the corresponding subplant curve due to the economic costassociated with resource consumption. High level optimization module 130can therefore be expected to select resource consumption values that lieon the corresponding subplant curve rather than above it.

The exception to this general rule is heat recovery chiller subplant 14.The equality constraints for heat recovery chiller subplant 14 providethat heat recovery chiller subplant 14 produces hot water at a rateequal to the subplant's cold water production plus the subplant'selectricity use. The inequality constraints generated by subplant curveincorporator 178 for heat recovery chiller subplant 14 allow high leveloptimization module 130 to overuse electricity to make more hot waterwithout increasing the amount of cold water production. This behavior isextremely inefficient and only becomes a realistic possibility when thedemand for hot water is high and cannot be met using more efficienttechniques. However, this is not how heat recovery chiller subplant 14actually operates.

To prevent high level optimization module 130 from overusingelectricity, subplant curve incorporator 178 may check whether thecalculated amount of electricity use (determined by the optimizationalgorithm) for heat recovery chiller subplant 14 is above thecorresponding subplant curve. In some embodiments, the check isperformed after each iteration of the optimization algorithm. If thecalculated amount of electricity use for heat recovery chiller subplant14 is above the subplant curve, subplant curve incorporator 178 maydetermine that high level optimization module 130 is overusingelectricity. In response to a determination that high level optimizationmodule 130 is overusing electricity, subplant curve incorporator 178 mayconstrain the production of heat recovery chiller subplant 14 at itscurrent value and constrain the electricity use of subplant 14 to thecorresponding value on the subplant curve. High level optimizationmodule 130 may then rerun the optimization with the new equalityconstraints.

Still referring to FIG. 4 , high level optimization module 130 is shownto include a ground loop module 152 and a heat exchanger module 154. Insome embodiments, central plant 10 includes a heat exchanger configuredto transfer heat between hot water loop 24 and condenser water loop 28.In some embodiments, central plant 10 includes a ground loop that servesas heat rejection for chiller subplant 16 and/or heat extraction forheat recovery chiller subplant 14. Ground loop module 152 and heatexchanger module 154 may be configured to modify the optimizationproblem to account for heat transfer resulting from operation of theheat exchanger and/or the ground loop.

Ground loop module 152 may incorporate heat rejection to the ground loopinto the optimization problem by changing the amount of electricity andwater usage by chiller subplant 16. For example, for loadings up to theheat rejection capacity of the ground loop, chiller subplant 16 may usean additional amount of electricity to run the ground loop pumps. Theadditional electricity usage may be constant or may vary per unit offlow through the ground loop. The amount of water production of chillersubplant 16 may be constant regardless of whether the ground loop isused.

Ground loop module 152 and heat exchanger module 154 may incorporateheat extraction from the ground loop and heat transfer between hot waterloop 24 and condenser water loop 28 into the optimization problem in asimilar manner. For example, ground loop module 152 and heat exchangermodule 154 may use heat extraction from the ground loop and heattransfer between loops 24 and 28 to modify the load seen by the centralplant equipment. Ground loop module 152 may use the ground loop tocreate what appears as a false building load to the equipment, therebyallowing heat recovery chiller subplant 14 to operate as heat pumps whenthe building load does not add enough heat to the system. This outcomemay be optimal when the ratio between electricity prices and gas pricesis low such that it is less expensive to operate the ground loop and theheat exchanger using electricity than it would be to use natural gas togenerate heat in heater subplant 12.

Heat exchanger module 154 may use the heat exchanger to create whatappears to be a false hot water building load, thereby allowing heatrecovery chiller subplant 14 to operate as conventional chillers. Theexcess heat from heat recovery chiller subplant 14 may be transferredthrough the heat exchanger to condenser loop 28 and ultimately into theatmosphere or into the ground. In some embodiments, heat exchangermodule 154 operates the heat exchanger to prevent condenser loop frombecoming overloaded. For example, heat exchanger module 154 may limitthe total heat rejected to the capacity of condenser loop 28 minus theheat produced by the conventional chillers.

Ground loop module 152 and heat exchanger module 154 may modify thedecision matrix x by adding a new decision vector for each type ofthermal energy load. The new decision vectors may represent theoverproduction of each thermal energy load for each time step k withinthe optimization period. For example, the modified decision matrix mayappear as follows:x=[{dot over (Q)} _(Chiller,1 . . . n) ,{dot over (Q)}_(hrChiller,1 . . . n) ,{dot over (Q)} _(Heater,1 . . . n) ,{dot over(Q)} _(HotStorage,1 . . . n) ,{dot over (Q)} _(ColdStorage,1 . . . n) ,. . . , Q _(ColdUnmet,1 . . . n) ,Q _(HotUnmet,1 . . . n) {dot over (Q)}_(ColdOver,1 . . . n) {dot over (Q)} _(HotOver,1 . . . n)]^(T)where {dot over (Q)}_(ColdOver,1 . . . n) and {dot over(Q)}_(HotOver,1 . . . n) are n-dimensional vectors representing theoverproduction rates of the cold thermal energy load and the hot thermalenergy load, respectively, for each time step k within the optimizationperiod.

Ground loop module 152 and heat exchanger module 154 may modify theequality constraints to account for any overproduced thermal energyloads. The overproduced thermal energy loads may be added to theequality constraints as slack variables that operate in the oppositedirection of the unmet loads. The modified equality constraints mayrequire that the predicted thermal energy loads plus any overproductionare equal to the total loads satisfied by subplants 12-22 plus anyunsatisfied thermal energy loads. The modified equality constraints canbe placed in the form Hx=g by defining the H matrix and the g vector asfollows:

${A = \begin{bmatrix}\left\lbrack I_{h} \right\rbrack & \left\lbrack I_{h} \right\rbrack & \left\lbrack 0_{h} \right\rbrack & \left\lbrack 0_{h} \right\rbrack & \left\lbrack I_{h} \right\rbrack & {\left\lbrack I_{h} \right\rbrack - \left\lbrack D_{- 1} \right\rbrack} & \left\lbrack 0_{h} \right\rbrack & {- \left\lbrack I_{h} \right\rbrack} & \left\lbrack 0_{h} \right\rbrack \\\left\lbrack 0_{h} \right\rbrack & \left( {1 + {u_{{elec},{hrChiller}}\left\lbrack I_{h} \right\rbrack}} \right. & \left\lbrack I_{h} \right\rbrack & \left\lbrack I_{h} \right\rbrack & \left\lbrack 0_{h} \right\rbrack & \left\lbrack 0_{h} \right\rbrack & {\left\lbrack I_{h} \right\rbrack - \left\lbrack D_{- 1} \right\rbrack} & \left\lbrack 0_{h} \right\rbrack & {- \left\lbrack I_{h} \right\rbrack}\end{bmatrix}}{g = \begin{bmatrix}{\hat{\ell}}_{{Cold},{1\ldots k}} \\{\hat{\ell}}_{{Hot},{1\ldots k}}\end{bmatrix}}$where [D⁻¹] is a lower diagonal matrix of ones. Ground loop module 152and heat exchanger module 154 may modify the cost vector c with theadditional cost of the pumping power per unit of overproduction requiredto run the ground loop and/or the heat exchanger.

Still referring to FIG. 4 , high level optimization module 130 is shownto include a demand charge module 156. As discussed above, optimizationframework module 142 may formulate the optimization problem as:

${{\underset{x}{\arg\min}c^{T}x};{{{subject}{to}{Ax}} \leq b}},{{Hx} = g}$However, such a formulation does not account for the demand charge.

The demand charge is an additional charge imposed by some utilityproviders based on the maximum rate of energy consumption during anapplicable demand charge period. For example, the demand charge may beprovided in terms of dollars per unit of power (e.g., $/kW) and may bemultiplied by the peak power usage (e.g., kW) during a demand chargeperiod to calculate the demand charge. In some instances, the demandcharge can account for more than 15% of the electrical bill. Failure toinclude the demand charge in the optimization scheme can cause all ofthe equipment to turn on at the same time (e.g., the most efficient orlowest cost time). This would be optimal from a consumption coststandpoint. However, shifting some of the load in time may savethousands of dollars on demand while only costing a few dollars inconsumption cost.

Demand charge module 156 may be configured to modify the optimizationproblem to account for the demand charge. Incorporating the demandcharge into the optimization framework may greatly improve theperformance of the high level optimization. For example, including thedemand charge in the optimization framework may reduce the totaloperating costs of central plant 10 by an additional 5% on top of the8-10% cost reduction provided by other modules of central plantcontroller 102. In various implementations, the savings provided bydemand charge module 156 and/or central plant controller 102 as a wholemay be greater than or less than the exemplary amounts defined hereindue to differences in plant configuration and/or energy costs.

Demand charge module 156 may account for the demand charge by modifyingthe cost function used by high level optimization module 130. Themodified cost function may be defined as:

${{\underset{x}{\arg\min}\left\lbrack {{c^{T}x} + {c_{demand}{\max\left( P_{{elec},k} \right)}}} \right\rbrack};{{{subject}{to}{Ax}} \leq b}},{{Hx} = g}$where c_(demand) is the demand charge (e.g., $/kW) for the applicabledemand charge period and P_(elec,k) is the total electrical powerconsumption of central plant 10 and the building/campus at time step k.The term max(P_(elec,k)) selects the peak electrical power consumptionat any time during the demand charge period. The demand chargec_(demand) and the demand charge period may be defined by the utilityrate information received from utilities 126 and may be provided to highlevel optimization module 130 by load/rate prediction module 122.

Incorporating the demand charge into the optimization frameworkcomplicates the optimization problem in two primary ways. First, thecost function is no longer linear due to the inclusion of the max( )function. Second, the consumption term c^(T)x calculates cost over aconsumption period defined by a time horizon, whereas the demand chargeterm c_(demand) max(P_(elec,k)) calculates cost over the demand chargeperiod. For example, the consumption period may be defined as the timeperiod beginning at the current time step k and ending at a future timestep k+h, where h represents the time horizon. The demand charge periodmay be defined by utilities 126 and provided to high level optimizationmodule 130 along with the utility rate information. In some instances,the consumption period and the demand charge period may not be the same.This complicates the optimization problem by obfuscating potentialtrade-offs between control decisions that reduce the consumption term atthe expense of the demand charge term or vice versa.

Demand charge module 156 may modify the optimization problem toincorporate the demand charge term into the linear optimizationframework. For example, demand charge module 156 may modify the decisionmatrix x by adding a new decision variable x_(peak) as follows:x _(new)=[ . . . u _(Chiller,elec,1 . . . n) . . . u_(hpChiller,elec,1 . . . n) . . . u _(Heater,elec,1 . . . n) . . . x_(peak)]^(T)where x_(peak) is the peak power consumption within the demand chargeperiod. Demand charge module 156 may modify the cost vector c asfollows:c _(new)=[ . . . c _(elec,1 . . . n) . . . c _(elec,1 . . . n) . . . c_(elec,1 . . . n) . . . c _(demand)]^(T)such that the demand charge c_(demand) is multiplied by the peak powerconsumption x_(peak).

Demand charge module 156 may formulate and/or apply inequalityconstraints to ensure that the peak power consumption x_(peak) isgreater than or equal to the maximum electric demand over the demandcharge period. I.e.:x _(peak)≥max(u _(Chiller,elec,k) +u _(hpChiller,elec,k) +u_(Heater,elec,k) +P _(elec,campus,k))∀k∈horizonThis inequality constraint may be represented in the linear optimizationframework by defining the A matrix and the b vector as follows:A=[ . . . [I _(h)] . . . [I _(h)] . . . [I _(h)] . . . −1], b=−P_(elec,campus,k)During the high level optimization process, high level optimizationmodule 130 may choose a x_(peak) that is equal to the maximum electricaldemand over the demand charge period to minimize the cost associatedwith x_(peak).

Demand charge module 156 may apply an inequality constraint to ensurethat the peak power consumption decision variable x_(peak) is greaterthan or equal to its previous value x_(peak),previous during the demandcharge period. This inequality constraint may be represented in thelinear optimization framework by defining the A matrix and the b vectoras follows:A=[ . . . −1], b=−x _(peak,previous)

Advantageously, the modifications to the decision variable matrix x, thecost vector c, and the inequality constraints provided by demand chargemodule 156 allow the cost function to be written in a linear form asfollows:

${{{\underset{x}{\arg\min}\left\lbrack {c_{new}^{T}x_{new}} \right\rbrack} = {\underset{x}{\arg\min}\left\lbrack {{c^{T}x} + {c_{demand}x_{peak}}} \right\rbrack}};{{{subject}{to}{Ax}} \leq b}},{{Hx} = g}$This linear form of the cost function can be used in the linearoptimization framework.

The cost function as written in the previous equation has componentsthat are over different time periods. For example, the consumption termc^(T)x is over the consumption period whereas the demand charge termc_(demand)x_(peak) is over the demand charge period. To properly makethe trade-off between increasing the demand charge versus increasing thecost of energy consumption, demand charge module 156 may apply aweighting factor to the demand charge term and/or the consumption term.For example, demand charge module 156 may divide the consumption termc^(T)x by the duration h of the consumption period (i.e., the timeperiod between the current time and the time horizon) and multiply bythe amount of time d_(demand) remaining in the current demand chargeperiod so that the entire cost function is over the demand chargeperiod. The new optimization function may be given by:

${{\underset{x}{\arg\min}\left\lbrack {{\frac{d_{demand}}{h}c^{T}x} + {c_{demand}x_{peak}}} \right\rbrack};{{{subject}{to}{Ax}} \leq b}},{{Hx} = g}$which is equivalent to:

${{\underset{x}{\arg\min}\left\lbrack {{c^{T}x} + {\frac{h}{d_{demand}}c_{demand}x_{peak}}} \right\rbrack};{{{subject}{to}{Ax}} \leq b}},{{Hx} = g}$The latter form of the new optimization function has the advantage ofadjusting only one term of the function rather than several.

Still referring to FIG. 4 , high level optimization module 130 is shownto include a load change penalty module 158. In some instances, highlevel optimization module 130 determines a solution to the optimizationproblem that includes significantly changing the load on one or more ofsubplants 12-22 within a relatively short timeframe. For example, thelowest cost solution from a resource consumption standpoint may involvetaking a subplant from off to full load and back to off again withinonly a few time steps. This behavior may result from high leveloptimization module 130 identifying small fluctuations in the economiccost of resources and operating central plant 10 accordingly to achievethe minimal economic cost. However, operating central plant 10 in such away may be undesirable due to various negative effects of rapidlychanging the subplant loads (e.g., increased equipment degradation),especially if the cost saved is relatively minimal (e.g., a few cents ordollars).

Load change penalty module 158 may modify the optimization problem tointroduce a penalty for rapidly changing the subplant loads. Forexample, load change penalty module 158 may modify the decision matrix xby adding a new decision vector for each subplant. The new decisionvectors represent the change in subplant load for each subplant from onetime step to the next. For example, load change penalty module 158 maymodify the decision matrix x as follows:x=[ . . . {dot over (Q)} _(Chiller,1 . . . n) . . . {dot over (Q)}_(hrChiller,1 . . . n) . . . {dot over (Q)} _(Heater,1 . . . n) . . .δ_(Chiller,1 . . . n) δ_(hrChiller,1 . . . n) δ_(Heater,1 . . . n)]^(T)where δ_(Chiller,1 . . . n), δ_(hrChiller,1 . . . n), andδ_(Heater,1 . . . n) are n-dimensional vectors representing the changein subplant load for {dot over (Q)}_(Chiller,1 . . . n), {dot over(Q)}_(hrChiller,1 . . . n), and {dot over (Q)}_(Heater,1 . . . n),respectively, at each time step k relative to the previous time stepk−1.

Load change penalty module 158 may modify the cost vector c to add acost associated with changing the subplant loads. For example, loadchange penalty module 158 may modify the cost vector c as follows:c=[ . . . 0_(n) . . . 0_(n) . . . 0_(n) . . . c _(δChiller,1 . . . n) c_(δhrChiller,1 . . . n) c _(δHeater,1 . . . n)]^(T)

Load change penalty module 158 may add constraints such that each of theload change variables δ cannot be less than the change in thecorresponding subplant load {dot over (Q)}. For example, the addedconstraints for chiller subplant 16 may have the following form:

${A = \begin{bmatrix}\cdots & {I_{h} - D_{- 1}} & \cdots & {- \left\lbrack I_{h} \right\rbrack} & \left\lbrack 0_{h} \right\rbrack & \left\lbrack 0_{h} \right\rbrack \\\cdots & {D_{- 1} - I_{h}} & \cdots & {- \left\lbrack I_{h} \right\rbrack} & \left\lbrack 0_{h} \right\rbrack & \left\lbrack 0_{h} \right\rbrack \\ \vdots & \vdots & \vdots & \vdots & \vdots & \vdots \end{bmatrix}},{b = \begin{bmatrix}\begin{bmatrix}{\overset{.}{Q}}_{{Chiller},{old}} \\0_{h - 1}\end{bmatrix} \\\begin{bmatrix}{- {\overset{.}{Q}}_{{Chiller},{old}}} \\0_{h - 1}\end{bmatrix} \\ \vdots \end{bmatrix}}$where {dot over (Q)}_(Chiller,old) is the value for {dot over(Q)}_(Chiller) at the previous time step. Similar constraints may beadded for each of subplants 12-22.

The constraints added by load change penalty module 158 require that theload change variables δ are greater than or equal to the magnitude ofthe difference between the current value of the corresponding subplantload {dot over (Q)} and the previous value of the subplant load {dotover (Q)}_(old). In operation, high level optimization module 130 mayselect values for the load change variables δ that are equal to themagnitude of the difference due to the costs associated with the loadchange variables. In other words, high level optimization module 130 maynot choose to make the load change variables δ greater than the actualchange in the corresponding subplant load because making the load changevariables δ greater than necessary would be suboptimal.

Still referring to FIG. 4 , high level optimization module 130 is shownto include a tank forced full module 160. Tank forced full module 160may modify the optimization problem such that the thermal energy storage(TES) tanks are forced to full at the end of the optimization period.This feature provides increased robustness in the event of a subplantfailure and/or controller failure by ensuring that the TES tanks havesufficient stored thermal energy to satisfy building loads while thefailure is being repaired. For example, plant operators can use thestored thermal energy to meet the building loads while central plantcontroller 102 is brought back online.

Tank forced full module 160 may force the TES tanks to full byincreasing the cost of discharging the TES tanks. In some embodiments,tank forced full module 160 modifies the cost of discharging the TEStanks such that the discharge cost is higher than other costs in thecost function, but less than the cost of unmet loads. This forces highlevel optimization module 130 to take the benefit (i.e., negative cost)of charging the TES tanks to their maximum values.

Referring now to FIGS. 5A-B, two subplant curves 500 and 510 are shown,according to an exemplary embodiment. Subplant curves 500 and 510 may beused as subplant curves 140, as described with reference to FIG. 3 .Subplant curves 500 and 510 define the resource usage of a subplant(e.g., one of subplants 12-22) as a function of the subplant load. Eachsubplant curve may be specific to a particular subplant and a particulartype of resource used by the subplant. For example, subplant curve 500may define the electricity use 502 of chiller subplant 16 as a functionof the load 504 on chiller subplant 16, whereas subplant curve 510 maydefine the water use 506 of chiller subplant 16 as a function of theload 504 on chiller subplant 16. Each of subplants 12-22 may have one ormore subplant curves (e.g., one for each type of resource consumed bythe subplant).

In some embodiments, low level optimization module 132 generatessubplant curves 500 and 510 based on equipment models 120 (e.g., bycombining equipment models 120 for individual devices into an aggregatecurve for the subplant). Low level optimization module 132 may generatesubplant curves 500 and 510 by running the low level optimizationprocess for several different loads and weather conditions to generatemultiple data points. Low level optimization module 132 may fit a curveto the data points to generate the subplant curves. In otherembodiments, low level optimization module 132 provides the data pointsto high level optimization module 132 and high level optimization module132 generates the subplant curves using the data points.

Referring now to FIG. 6 , another subplant curve 600 is shown, accordingto an exemplary embodiment. Subplant curve 600 defines the electricityuse of chiller subplant 16 (i.e., u_(Chiller,elec)) as a function of thecold water production of chiller subplant 16 (i.e., {dot over(Q)}_(Chiller)). In some embodiments, subplant curve 600 is generated bycombining efficiency curves for individual devices of chiller subplant16 (e.g., individual chillers, pumps, etc.). For example, each of thechillers in subplant 16 may have a device-specific efficiency curve thatdefines the amount of electricity use by the chiller as a function ofthe load on the chiller. Many devices operate less efficiently at higherloads and have device efficiency curves that are U-shaped functions ofload. Accordingly, combining multiple device efficiency curves to formsubplant curve 600 may result in subplant curve 600 having one or morewaves 602, as shown in FIG. 6 . Waves 602 may be caused by a singledevice loading up before it is more efficient to turn on another deviceto satisfy the subplant load.

Referring now to FIG. 7 , a linearized subplant curve 700 is shown,according to an exemplary embodiment. Subplant curve 700 defines theelectricity use of chiller subplant 16 (i.e., u_(Chiller,elec)) as afunction of the cold water production of chiller subplant 16 (i.e., {dotover (Q)}_(Chiller)). Subplant curve 700 may be generated by convertingsubplant curve 600 into a linearized convex curve. A convex curve is acurve for which a line connecting any two points on the curve is alwaysabove or along the curve (i.e., not below the curve). Convex curves maybe advantageous for use in the high level optimization because theyallow for an optimization process that is less computationally expensiverelative to an optimization process that uses non-convex functions.

In some embodiments, subplant curve 700 is generated by subplant curvelinearizer 176, as described with reference to FIG. 4 . Subplant curve700 may be created by generating a plurality of linear segments (i.e.,segments 702, 704, and 706) that approximate subplant curve 600 andcombining the linear segments into a piecewise-defined linearized convexcurve 700. Linearized subplant curve 700 is shown to include a firstlinear segment 702 connecting point [u₁, Q₁] to point [u₂, Q₂], a secondlinear segment 704 connecting point [u₂, Q₂] to point [u₃, Q₃], and athird linear segment 706 connecting point [u₃, Q₃] to point [u₄, Q₄].The endpoints of line segments 702-706 may be used to form constraintsthat force the electricity use of chiller subplant 16 in the epigraph ofthe linearized subplant curve 700.

Referring now to FIG. 8 , another subplant curve 800 is shown, accordingto an exemplary embodiment. Subplant curve 800 defines the energy use ofone of subplants 12-22 as a function of the load on the subplant forseveral different weather conditions. In some embodiments, subplantcurve 800 is generated by subplant curves module 170 using experimentaldata obtained from the low level optimization module 132. For example,subplant curve updater 172 may request resource usage data from lowlevel optimization module 132 for various combinations of loadconditions and environmental conditions. In the embodiment shown in FIG.8 , subplant curve updater 172 requests energy use data for eachcombination of temperature (e.g., 40° F., 50° F., 60° F., and 70° F.)and load (e.g., 170 tons, 330 tons, 500 tons, 830 tons, and 1000 tons).Low level optimization module 132 may perform the low level optimizationprocess for the requested load and temperature combinations and returnan energy use value for each combination.

Subplant curve updater 172 may use the data points provided by low leveloptimization module 132 to find the best piecewise linear convexfunction that fits the data. For example, subplant curve updater 172 mayfit a first subplant curve 802 to the data points at 70° F., a secondsubplant curve 804 to the data points at 60° F., a third subplant curve806 to the data points at 50° F., and a fourth subplant curve 808 to thedata points at 40° F. Subplant curve updater 172 may store the generatedsubplant curves 802-808 in subplant curves database 174 for use in thehigh level optimization algorithm.

In some implementations, central plant controller 102 uses high leveloptimization module 130 as part of an operational tool to exercisereal-time control over central plant 10. In the operational tool, highlevel optimization module 130 may receive load and rate predictions fromload/rate prediction module 122 and subplant curves (or data that can beused to generate subplant curves) from low level optimization module132. When implemented in the operational tool, high level optimizationmodule 130 may determine optimal subplant loads for heater subplant 12,heat recovery chiller subplant 14, chiller subplant 16, hot TES subplant20, and/or cold TES subplant 22, as described with reference to FIGS.2-4 . In some embodiments, high level optimization module 130 determinesground loop and heat exchanger transfer rates in addition to thesubplant loads. When implemented in the operational tool, high leveloptimization module 130 may provide the determined subplant loads and/ortransfer rates to low level optimization module 132 for use indetermining optimal on/off decisions and/or operating setpoints for theequipment of each subplant.

Referring now to FIG. 9 , a block diagram of a planning system 900 isshown, according to an exemplary embodiment. Planning system 900 may beconfigured to use optimization module 128 as part of a planning tool 902to simulate the operation of a central plant over a predetermined timeperiod (e.g., a day, a month, a week, a year, etc.) for planning,budgeting, and/or design considerations. When implemented in planningtool 902, optimization module 128 may operate in a similar manner asdescribed with reference to FIGS. 2-4 . For example, optimization module128 may use building loads and utility rates to determine an optimalsubplant load distribution to minimize cost over a simulation period.However, planning tool 902 may not be responsible for real-time controlof a building automation system or central plant.

In planning tool 902, high level optimization module 130 may receiveplanned loads and utility rates for the entire simulation period. Theplanned loads and utility rates may be defined by input received from auser via a client device 922 (e.g., user-defined, user selected, etc.)and/or retrieved from a plan information database 926. High leveloptimization module 130 uses the planned loads and utility rates inconjunction with subplant curves from low level optimization module 132to determine optimal subplant loads (i.e., an optimal dispatch schedule)for a portion of the simulation period.

The portion of the simulation period over which high level optimizationmodule 130 optimizes the subplant loads may be defined by a predictionwindow ending at a time horizon. With each iteration of theoptimization, the prediction window is shifted forward and the portionof the dispatch schedule no longer in the prediction window is accepted(e.g., stored or output as results of the simulation). Load and ratepredictions may be predefined for the entire simulation and may not besubject to adjustments in each iteration. However, shifting theprediction window forward in time may introduce additional planinformation (e.g., planned loads and/or utility rates) for thenewly-added time slice at the end of the prediction window. The new planinformation may not have a significant effect on the optimal dispatchschedule since only a small portion of the prediction window changeswith each iteration.

In some embodiments, high level optimization module 130 requests all ofthe subplant curves used in the simulation from low level optimizationmodule 132 at the beginning of the simulation. Since the planned loadsand environmental conditions are known for the entire simulation period,high level optimization module 130 may retrieve all of the relevantsubplant curves at the beginning of the simulation. In some embodiments,low level optimization module 132 generates functions that map subplantproduction to equipment level production and resource use when thesubplant curves are provided to high level optimization module 130.These subplant to equipment functions may be used to calculate theindividual equipment production and resource use (e.g., in apost-processing module) based on the results of the simulation.

Still referring to FIG. 9 , planning tool 902 is shown to include acommunications interface 904 and a processing circuit 906.Communications interface 904 may include wired or wireless interfaces(e.g., jacks, antennas, transmitters, receivers, transceivers, wireterminals, etc.) for conducting data communications with varioussystems, devices, or networks. For example, communications interface 904may include an Ethernet card and port for sending and receiving data viaan Ethernet-based communications network and/or a WiFi transceiver forcommunicating via a wireless communications network. Communicationsinterface 904 may be configured to communicate via local area networksor wide area networks (e.g., the Internet, a building WAN, etc.) and mayuse a variety of communications protocols (e.g., BACnet, IP, LON, etc.).

Communications interface 904 may be a network interface configured tofacilitate electronic data communications between planning tool 902 andvarious external systems or devices (e.g., client device 922, resultsdatabase 928, plan information database 926, etc.). For example,planning tool 902 may receive planned loads and utility rates fromclient device 922 and/or plan information database 926 viacommunications interface 904. Planning tool 902 may use communicationsinterface 904 to output results of the simulation to client device 922and/or to store the results in results database 928.

Still referring to FIG. 9 , processing circuit 906 is shown to include aprocessor 910 and memory 912. Processor 910 may be a general purpose orspecific purpose processor, an application specific integrated circuit(ASIC), one or more field programmable gate arrays (FPGAs), a group ofprocessing components, or other suitable processing components.Processor 910 may be configured to execute computer code or instructionsstored in memory 912 or received from other computer readable media(e.g., CDROM, network storage, a remote server, etc.).

Memory 912 may include one or more devices (e.g., memory units, memorydevices, storage devices, etc.) for storing data and/or computer codefor completing and/or facilitating the various processes described inthe present disclosure. Memory 912 may include random access memory(RAM), read-only memory (ROM), hard drive storage, temporary storage,non-volatile memory, flash memory, optical memory, or any other suitablememory for storing software objects and/or computer instructions. Memory912 may include database components, object code components, scriptcomponents, or any other type of information structure for supportingthe various activities and information structures described in thepresent disclosure. Memory 912 may be communicably connected toprocessor 910 via processing circuit 906 and may include computer codefor executing (e.g., by processor 906) one or more processes describedherein.

Still referring to FIG. 9 , memory 912 is shown to include a GUI engine926, web services 914, and configuration tools 918. In an exemplaryembodiment, GUI engine 916 includes a graphical user interface componentconfigured to provide graphical user interfaces to a user for selectingor defining plan information for the simulation (e.g., planned loads,utility rates, environmental conditions, etc.). Web services 914 mayallow a user to interact with planning tool 902 via a web portal and/orfrom a remote system or device (e.g., an enterprise controlapplication).

Configuration tools 918 can allow a user to define (e.g., via graphicaluser interfaces, via prompt-driven “wizards,” etc.) various parametersof the simulation such as the number and type of subplants, the deviceswithin each subplant, the subplant curves, device-specific efficiencycurves, the duration of the simulation, the duration of the predictionwindow, the duration of each time step, and/or various other types ofplan information related to the simulation. Configuration tools 918 canpresent user interfaces for building the simulation. The user interfacesmay allow users to define simulation parameters graphically. In someembodiments, the user interfaces allow a user to select a pre-stored orpre-constructed simulated plant and/or plan information (e.g., from planinformation database 926) and adapt it or enable it for use in thesimulation.

Still referring to FIG. 9 , memory 912 is shown to include optimizationmodule 128. Optimization module 128 may use the planned loads andutility rates to determine optimal subplant loads over a predictionwindow. The operation of optimization module 128 may be the same orsimilar as previously described with reference to FIGS. 2-4 . With eachiteration of the optimization process, optimization module 128 may shiftthe prediction window forward and apply the optimal subplant loads forthe portion of the simulation period no longer in the prediction window.Optimization module 128 may use the new plan information at the end ofthe prediction window to perform the next iteration of the optimizationprocess. Optimization module 128 may output the applied subplant loadsto reporting applications 930 for presentation to a client device 922(e.g., via user interface 924) or storage in results database 928.

Still referring to FIG. 9 , memory 912 is shown to include reportingapplications 930. Reporting applications 930 may receive the optimizedsubplant loads from optimization module 128 and, in some embodiments,costs associated with the optimized subplant loads. Reportingapplications 930 may include a web-based reporting application withseveral graphical user interface (GUI) elements (e.g., widgets,dashboard controls, windows, etc.) for displaying key performanceindicators (KPI) or other information to users of a GUI. In addition,the GUI elements may summarize relative energy use and intensity acrossvarious plants, subplants, or the like. Other GUI elements or reportsmay be generated and shown based on available data that allow users toassess the results of the simulation. The user interface or report (orunderlying data engine) may be configured to aggregate and categorizesubplant loads and the costs associated therewith and provide theresults to a user via a GUI. The GUI elements may include charts orhistograms that allow the user to visually analyze the results of thesimulation. An exemplary output that may be generated by reportingapplications 930 is shown in FIG. 10 .

Referring now to FIG. 10 , several graphs 1000 illustrating theoperation of planning tool 902 are shown, according to an exemplaryembodiment. With each iteration of the optimization process, planningtool 902 selects an optimization period (i.e., a portion of thesimulation period) over which the optimization is performed. Forexample, planning tool 902 may select optimization period 1002 for usein the first iteration. Once the optimal load distribution 1010 has beendetermined, planning tool 902 may select a portion 1018 of loaddistribution 1010 to send to plant dispatch 1030. Portion 1018 may bethe first b time steps of load distribution 1010. Planning tool 902 mayshift the optimization period 1002 forward in time, resulting inoptimization period 1004. The amount by which the prediction window isshifted may correspond to the duration of time steps b.

Planning tool 902 may repeat the optimization process for optimizationperiod 1004 to determine the optimal plant load distribution 1012.Planning tool 902 may select a portion 1020 of plant load distribution1012 to send to plant dispatch 1030. Portion 1020 may be the first btime steps of load distribution 1012. Planning tool 902 may then shiftthe prediction window forward in time, resulting in optimization period1006. This process may be repeated for each subsequent optimizationperiod (e.g., optimization periods 1006, 1008, etc.) to generate updatedload distributions (e.g., load distributions 1014, 1016, etc.) and toselect portions of each load distribution (e.g., portions 1022, 1024) tosend to plant dispatch 1030. Plant dispatch 1030 includes the first btime steps 1018-1024 from each of optimization periods 1002-1008. Oncethe optimal subplant load distribution 1030 is compiled for the entiresimulation period, the results may be sent to reporting applications930, results database 928, and/or client device 922, as described withreference to FIG. 9 .

Referring now to FIG. 11 , a flowchart of a process 1100 for optimizingcost in a central plant is shown, according to an exemplary embodiment.In various implementations, process 1100 may be performed by centralplant controller 102 or planning tool 902. The central plant may includea plurality of subplants (e.g., subplants 12-22) configured to serve theenergy loads of a building or campus. The central plant may be an actualplant (e.g., central plant 10) or a simulated central plant including aplurality of simulated subplants.

Process 1100 is shown to include receiving load prediction data andutility rate data (step 1102). The load prediction data may includepredicted or planned thermal energy loads for a building or campus foreach time step k (e.g., k=1 . . . n) of an optimization period. The loadprediction data may include predicted or planned values one or moredifferent types of loads for the building or campus. For example, theload prediction data may include a predicted hot water load

and a predicted cold water load

for each time step k within the prediction window.

In some embodiments, the load prediction data are based on weatherforecasts from a weather service and/or feedback from the building orcampus. Feedback from the building or campus may include various typesof sensory inputs (e.g., temperature, flow, humidity, enthalpy, etc.) orother data relating to the controlled building (e.g., inputs from a HVACsystem, a lighting control system, a security system, a water system,etc.). In some embodiments, the load prediction data are generated byload/rate prediction module 122, as described with reference to FIG. 2 .For example, the load prediction data may be based on a measuredelectric load and/or previous measured load data from the building orcampus. The load prediction data may be a function of a given weatherforecast ({circumflex over (ϕ)}_(w)), a day type (day), the time of day(t), and/or previous measured load data (Y_(k-1)). Such a relationshipis expressed in the following equation:

=f(ϕ_(w),day,t|Y _(k-1))

The utility rate data may indicate a cost or price per unit of one ormore resources (e.g., electricity, natural gas, water, etc.) consumed bythe central plant to serve the thermal energy loads of the building orcampus at each time step k in the prediction window. In someembodiments, the utility rates are time-variable rates. For example, theprice of electricity may be higher at certain times of day or days ofthe week (e.g., during high demand periods) and lower at other times ofday or days of the week (e.g., during low demand periods). The utilityrates may define various time periods and a cost per unit of a resourceduring each time period. Utility rates may be actual rates (e.g.,received from utilities 126) or predicted utility rates (e.g., estimatedby load/rate prediction module 122).

In some embodiments, the utility rates include demand charges for one ormore of the resources consumed by the central plant. A demand charge maydefine a separate cost based on the maximum usage of a particularresource (e.g., maximum energy consumption) during a demand chargeperiod. The utility rates may define various demand charge periods andone or more demand charges associated with each demand charge period. Insome instances, demand charge periods may overlap partially orcompletely with each other and/or with the prediction window. Theutility rate data may include time-variable (e.g., hourly) prices, amaximum service level (e.g., a maximum rate of consumption allowed bythe physical infrastructure or by contract) and, in the case ofelectricity, a demand charge or a charge for the peak rate ofconsumption within a certain period.

Still referring to FIG. 11 , process 1100 is shown to include generatingan objective function that expresses a total monetary cost of operatingthe central plant as a function of the utility rate data and an amountof resources consumed by the central plant (step 1104). In someembodiments, the objective function is a high level cost function J_(HL)for the central plant. The high level cost function J_(HL) may representthe sum of the monetary costs of each utility consumed by the centralplant for the duration of the optimization period. For example, the highlevel cost function J_(HL) may be described using the followingequation:

${J_{HL}\left( \theta_{HL} \right)} = {\sum\limits_{k = 1}^{n_{h}}{\sum\limits_{i = 1}^{n_{s}}\left\lbrack {\sum\limits_{j = 1}^{n_{u}}{{t_{s} \cdot c_{jk}}{u_{jik}\left( \theta_{HL} \right)}}} \right\rbrack}}$where n_(h) is the number of time steps k in the optimization period,n_(s) is the number of subplants, t_(s) is the duration of a time step,c_(jk) is the economic cost of utility j at a time step k of theoptimization period, and u_(jik) is the rate of use of utility j bysubplant i at time step k.

In some embodiments, the objective function is generated using a linearprogramming framework. For example, step 1104 may include generating anobjective function of the form:

${{\underset{x}{\arg\min}c^{T}x};{{{subject}{to}{Ax}} \leq b}},{{Hx} = g}$where c is a cost vector, x is a decision matrix, A and b are a matrixand vector (respectively) which describe inequality constraints on thevariables in the decision matrix x, and H and g are a matrix and vector(respectively) which describe equality constraints on the variables inthe decision matrix x. In other embodiments, the objective function maybe generated using any of a variety of other optimization frameworks(e.g., quadratic programming, linear-fractional programming, nonlinearprogramming, combinatorial algorithms, etc.).

In some embodiments, step 1104 includes formulating the decision matrixx. The loads across each of subplants 12-22 may be the decisionvariables in the decision matrix x. For example, for a central plantthat includes chillers, heat recovery chillers, hot water generators,and thermal energy storage, step 1104 may include formulating thedecision matrix x as:x=[{dot over (Q)} _(Chiller,1 . . . n) , {dot over (Q)}_(hrChiller,1 . . . n) , {dot over (Q)} _(Heater,1 . . . n) , {dot over(Q)} _(HotStorage,1 . . . n) , {dot over (Q)}_(ColdStorage,1 . . . n)]^(T)where {dot over (Q)}_(Chiller,1 . . . n), {dot over(Q)}_(hrChiller,1 . . . n), {dot over (Q)}_(Heater,1 . . . n), {dot over(Q)}_(HotStorage,1 . . . n), and {dot over (Q)}_(ColdStorage,1 . . . n)are n-dimensional vectors representing the thermal energy load assignedto chiller subplant 16, heat recovery chiller subplant 14, heatersubplant 12, hot TES subplant 20, and cold TES subplant 22,respectively, for each of the n time steps within the optimizationperiod.

In some embodiments, step 1104 includes generating the decision matrix xto include one or more decision vectors representing the resourceconsumption of each subplant. For example, for a central plant thatincludes a chiller subplant, step 1104 may include generating thedecision matrix x as follows:x=[ . . . {dot over (Q)} _(Chiller,1 . . . n) . . . u_(Chiller,elec,1 . . . n) u _(Chiller,water,1 . . . n) . . . ]^(T)where u_(Chiller,elec,1 . . . n) and u_(Chiller,water,1 . . . n) aren-dimensional vectors representing the amount of electrical consumptionand water consumption, respectively, by the chiller subplant at eachtime step k.

Step 1104 may include adding one or more resource consumption vectors tomatrix x for each of subplants 12-22. The decision vectors added in step1104 for a given subplant may represent an amount of resourceconsumption for each resource consumed by the subplant (e.g., water,electricity, natural gas, etc.) at each time step k within theoptimization period. For example, if a heater subplant consumes naturalgas, electricity, and water, step 1104 may include adding a decisionvector u_(Heater,gas,1 . . . n) representing an amount of natural gasconsumed by the heater subplant at each time step, a decision vectoru_(Heater,elec,1 . . . n) representing an amount of electricity consumedby the heater subplant at each time step, and a decision vectoru_(Heater,water,1 . . . n) representing an amount of water consumed bythe heater subplant at each time step. Step 1104 may include addingresource consumption vectors for other subplants in a similar manner.

In some embodiments, step 1104 includes generating the cost vector c.Generating the cost vector c may include adding economic costsassociated with the resource consumption required to produce thesubplant loads. For example, the decision matrix x provided above, step1104 may include generating the cost vector c as follows:c=[ . . . 0_(n) . . . c _(elec,1 . . . n) c _(water,1 . . . n) . . .]^(T)where 0_(n) is a n-dimensional zero vector indicating that the directeconomic cost of {dot over (Q)}_(Chiller,1 . . . n) is zero at each timestep, c_(elec,1 . . . n) is a n-dimensional vector indicating the perunit cost of electricity at each time step, and c_(water,1 . . . n) is an-dimensional vector indicating the per unit cost of water at each timestep. The cost vector associates an economic cost with the resourcesconsumed to produce the subplant loads rather than the subplant loadsthemselves. In some embodiments, the values for c_(elec,1 . . . n) andc_(water,1 . . . n) are utility rates obtained from the utility ratedata received in step 1102.

In some embodiments, step 1104 includes generating the A matrix and theb vector which describe the inequality constraints, and the H matrix andthe g vector which describe the equality constraints. The inequalityconstraints and equality constraints may be generated by inequalityconstraints module 146 and equality constraints module 148, as describedwith reference to FIG. 4 . For example, step 1104 may include generatinginequality constraints that constrain the decision variables in matrix xto be less than or equal to maximum capacities for the correspondingcentral plant equipment and less than or equal to maximumcharge/discharge rates for thermal energy storage. Step 1104 may includegenerating inequality constraints that prevent charging the thermalenergy storage above maximum capacity and/or discharging the thermalenergy storage below zero. Step 1104 may include generating equalityconstraints that ensure the building energy loads are satisfied at eachof the time steps in the prediction window.

In some embodiments, step 1104 includes modifying the objective functionto account for unmet loads (e.g., as described with reference to unmetloads module 150), to account for heat extraction or rejection to aground loop (e.g., as described with reference to ground loop module152), to account for heat exchange between the hot water loop and thecondenser water loop (e.g., as described with reference to heatexchanger module 154), to account for subplant curves that are notsimple linear functions of load (e.g., as described with reference tosubplant curves module 170), and/or to force the thermal energy storagetanks to full at the end of the prediction window (e.g., as describedwith reference to tank forced full module 160). Modifying the objectivefunction may include modifying the decision matrix x, the cost vector c,the A matrix and the b vector which describe the inequality constraints,and/or the H matrix and the g vector which describe the equalityconstraints.

Still referring to FIG. 11 , process 1100 is shown to include modifyingthe objective function to account for a demand charge (step 1106). Step1106 is an optional step that may be performed by demand charge module156 to account for a demand charge that may be imposed by utilityproviders in some pricing scenarios. The demand charge is an additionalcharge imposed by some utility providers based on the maximum rate ofenergy consumption during an applicable demand charge period. Forexample, the demand charge may be provided in terms of dollars per unitof power (e.g., $/kW) and may be multiplied by the peak power usage(e.g., kW) during a demand charge period to calculate the demand charge.

Accounting for the demand charge may include modifying the variouscomponents of the objective function such as the decision matrix x, thecost vector c, and/or the A matrix and the b vector which describe theinequality constraints. The modified objective function may be definedas:

${{\underset{x}{\arg\min}\left\lbrack {{c^{T}x} + {c_{demand}{\max\left( P_{{elec},k} \right)}}} \right\rbrack};{{{subject}{to}{Ax}} \leq b}},{{Hx} = g}$where c_(demand) is the demand charge for the applicable demand chargeperiod and P_(elec,k) is the total electrical power consumption of thecentral plant and the building/campus at time step k. The termmax(P_(elec,k)) selects the peak electrical power consumption at anytime during the demand charge period. The demand charge c_(demand) andthe demand charge period may be defined by the utility rate informationreceived in step 1102.

Step 1106 may include modifying the decision matrix x by adding a newdecision variable x_(peak) as follows:x _(new)=[ . . . u _(Chiller,elec,1 . . . n) . . . u_(hpChiller,elec,1 . . . n) . . . u _(Heater,elec,1 . . . n) . . . x_(peak)]^(T)where x_(peak) is the peak power consumption within the optimizationperiod. Step 1106 may include modifying the cost vector c as follows:c _(new)=[ . . . c _(elec,1 . . . n) . . . c _(elec,1 . . . n) . . . c_(elec,1 . . . n) . . . c _(demand)]^(T)such that the demand charge c_(demand) is multiplied by the peak powerconsumption x_(peak).

Step 1106 may include generating and/or imposing inequality constraintsto ensure that the peak power consumption x_(peak) is greater than orequal to the maximum electric demand for each time step in theoptimization period. I.e.:x _(peak)≥max(u _(Chiller,elec,k) +u _(hpChiller,elec,k) +u_(Heater,elec,k) +P _(elec,campus,k))∀k∈horizonThis inequality constraint may be represented in the linear optimizationframework by defining the A matrix and the b vector as follows:A=[ . . . [I _(h)] . . . [I _(h)] . . . [I _(h)] . . . −1], b=−P_(elec,campus,k)

Step 1106 may include generating and/or imposing an inequalityconstraint to ensure that the peak power consumption decision variablex_(peak) is greater than or equal to its previous valuex_(peak),previous during the demand charge period. This inequalityconstraint may be represented in the linear optimization framework bydefining the A matrix and the b vector as follows:A=[ . . . −1], b=−x _(peak,previous)

Advantageously, the modifications to the decision variable matrix x, thecost vector c, and the inequality constraints in step 1106 may allow theobjective function to be written in a linear form as follows:

${{{\underset{x}{\arg\min}\left\lbrack {c_{new}^{T}x_{new}} \right\rbrack} = {\underset{x}{\arg\min}\left\lbrack {{c^{T}x} + {c_{demand}x_{peak}}} \right\rbrack}};{{{subject}{to}{Ax}} \leq b}},{{Hx} = g}$This linear form of the objective function can be used in the linearoptimization framework.

In some embodiments, step 1106 includes applying a weighting factor toat least one of the consumption term and the demand charge term of theobjective function. For example, the objective function as written inthe previous equation has components that are over different timeperiods. The consumption term c^(T)x is over the consumption periodwhereas the demand charge term c_(demand)x_(peak) is over the demandcharge period. To properly make the trade-off between increasing thedemand charge versus increasing the cost of energy consumption, step1106 may include applying a weighting factor to the demand charge termand/or the consumption term. For example, step 1106 may include dividingthe consumption term c^(T)x by the duration h of the consumption period(i.e., the time period between the current time and the time horizon)and multiplying by the amount of time d_(demand) remaining in thecurrent demand charge period so that the entire objective function isover the demand charge period. The new optimization function may begiven by:

${{\underset{x}{\arg\min}\left\lbrack {{\frac{d_{demand}}{h}c^{T}x} + {c_{demand}x_{peak}}} \right\rbrack};{{{subject}{to}{Ax}} \leq b}},{{Hx} = g}$which is equivalent to:

${{\underset{x}{\arg\min}\left\lbrack {{c^{T}x} + {\frac{h}{d_{demand}}c_{demand}x_{peak}}} \right\rbrack};{{{subject}{to}{Ax}} \leq b}},{{Hx} = g}$The latter form of the new optimization function has the advantage ofadjusting only one term of the function rather than several.

Still referring to FIG. 11 , process 1100 is shown to include modifyingthe objective function to account for a load change penalty (step 1108).Step 1108 is an optional step that may be performed by load changepenalty module 158 to account for the cost of changing the loadsassigned to each of the subplants. In some instances, the lowest costsolution from a resource consumption standpoint may involve taking asubplant from off to full load and back to off again within only a fewtime steps. However, operating the central plant in such a way may beundesirable due to various negative effects of rapidly changing thesubplant loads (e.g., increased equipment degradation), especially ifthe cost saved is relatively minimal (e.g., a few cents or dollars).

Step 1108 may include modifying the objective function to introduce apenalty for rapidly changing the subplant loads. In some embodiments,step 1108 includes modifying the decision matrix x by adding a newdecision vector for each subplant. The new decision vectors representthe change in subplant load for each subplant from one time step to thenext. For example, step 1108 may include modifying the decision matrix xas follows:x=[ . . . {dot over (Q)} _(Chiller,1 . . . n) . . . {dot over (Q)}_(hrChiller,1 . . . n) . . . {dot over (Q)} _(Heater,1 . . . n) . . .δ_(Chiller,1 . . . n) δ_(hrChiller,1 . . . n) δ_(Heater,1 . . . n)]^(T)where δ_(Chiller,1 . . . n), δ_(hrChiller,1 . . . n), andδ_(Heater,1 . . . n) are n-dimensional vectors representing the changein subplant load for {dot over (Q)}_(Chiller,1 . . . n), {dot over(Q)}_(hrChiller,1 . . . n), and {dot over (Q)}_(Heater,1 . . . n),respectively, at each time step k relative to the previous time stepk-1.

Step 1108 may include modifying the cost vector c to add a costassociated with changing the subplant loads. For example, step 1108 mayinclude modifying the cost vector c as follows:c=[ . . . 0_(n) . . . 0_(n) . . . 0_(n) . . . c _(δChiller,1 . . . n) c_(δhrChiller,1 . . . n) c _(δHeater,1 . . . n)]^(T)

Step 1108 may include adding constraints such that each of the loadchange variables δ cannot be less than the change in the correspondingsubplant load {dot over (Q)}. For example, the added constraints for achiller subplant may have the following form:

${A = \begin{bmatrix}\cdots & {I_{h} - D_{- 1}} & \cdots & {- \left\lbrack I_{h} \right\rbrack} & \left\lbrack 0_{h} \right\rbrack & \left\lbrack 0_{h} \right\rbrack \\\cdots & {D_{- 1} - I_{h}} & \cdots & {- \left\lbrack I_{h} \right\rbrack} & \left\lbrack 0_{h} \right\rbrack & \left\lbrack 0_{h} \right\rbrack \\ \vdots & \vdots & \vdots & \vdots & \vdots & \vdots \end{bmatrix}},{b = \begin{bmatrix}\begin{bmatrix}{\overset{.}{Q}}_{{Chiller},{old}} \\0_{h - 1}\end{bmatrix} \\\begin{bmatrix}{- {\overset{.}{Q}}_{{Chiller},{old}}} \\0_{h - 1}\end{bmatrix} \\ \vdots \end{bmatrix}}$where {dot over (Q)}_(Chiller,old) is the value for {dot over(Q)}_(Chiller) at the previous time step. Similar constraints may beadded for each of subplants 12-22. The constraints added by in step 1108may require that the load change variables δ are greater than or equalto the magnitude of the difference between the current value of thecorresponding subplant load {dot over (Q)} and the previous value of thesubplant load {dot over (Q)}_(old).

Still referring to FIG. 11 , process 1100 is shown to include optimizingthe objective function over an optimization period subject to a set ofconstraints to determine an optimal distribution of energy loads overmultiple groups of central plant equipment (step 1110). The set ofconstraints may include the inequality constraints and the equalityconstraints formulated in steps 1104, 1106, and/or 1108. Optimizing theobjective function may include determining an optimal decision matrix x*that minimizes the cost function c^(T)x. The optimal decision matrix x*may correspond to the optimal decisions θ_(HL)* (for each time step kwithin an optimization period) that minimize the high level costfunction J_(HL), as described with reference to FIG. 3 .

Step 1110 may include using any of a variety of linear optimizationtechniques to determine the optimal decision matrix. For example, step1110 may include using basis exchange algorithms (e.g., simplex,crisscross, etc.), interior point algorithms (e.g., ellipsoid,projective, path-following, etc.), covering and packing algorithms,integer programming algorithms (e.g., cutting-plant, branch and bound,branch and cut, branch and price, etc.), or any other type of linearoptimization algorithm or technique to solve for the optimal decisionmatrix subject to the optimization constraints. For embodiments in whichnonlinear optimization is used, step 1110 may include using any of avariety of nonlinear optimization techniques to solve for the optimaldecision matrix. The result of step 1110 may be an optimal distributionof energy loads over the multiple groups of subplant equipment (i.e.,the multiple subplants) for each of the time steps k.

Still referring to FIG. 11 , process 1100 is shown to include using theoptimal distribution of energy loads to determine optimal operatingstatuses for individual devices of the central plant equipment (step1112). In some embodiments, step 1112 is performed by low leveloptimization module 132, as described with reference to FIGS. 2-4 . Forexample, step 1112 may include using the subplant loads determined instep 1110 to determine optimal low level decisions θ_(LL)* (e.g. binaryon/off decisions, flow setpoints, temperature setpoints, etc.) for thecentral plant equipment. In some embodiments, step 1112 is performed foreach of the plurality of subplants.

Step 1112 may include determining which devices of each subplant to useand/or the operating setpoints for such devices that will achieve thesubplant load setpoint while minimizing energy consumption. The lowlevel optimization performed in step 1112 may be described using thefollowing equation:

$\theta_{LL}^{*} = {\arg\min\limits_{\theta_{LL}}{J_{LL}\left( \theta_{LL} \right)}}$where θ_(LL)* contains the optimal low level decisions and LL is the lowlevel cost function.

To find the optimal low level decisions L, step 1112 may includeminimizing the low level cost function J_(LL). The low level costfunction J_(LL) may represent the total energy consumption for all ofthe equipment in the applicable subplant. The low level cost functionJ_(LL) may be described using the following equation:

${J_{LL}\left( \theta_{LL} \right)} = {\sum\limits_{j = 1}^{N}{t_{s} \cdot b_{j} \cdot {u_{j}\left( \theta_{LL} \right)}}}$where N is the number of devices of in the subplant, t_(s) is theduration of a time step, b_(j) is a binary on/off decision (e.g., 0=off,1=on), and u_(j) is the energy used by device j as a function of thesetpoint θ_(LL). Each device may have continuous variables which can bechanged to determine the lowest possible energy consumption for theoverall input conditions.

In some embodiments, step 1112 includes minimizing the low level costfunction J_(LL) subject to inequality constraints based on thecapacities of the subplant equipment and equality constraints based onenergy and mass balances. In some embodiments, the optimal low leveldecisions θ_(LL)* are constrained by switching constraints defining ashort horizon for maintaining a device in an on or off state after abinary on/off switch. The switching constraints may prevent devices frombeing rapidly cycled on and off.

Step 1112 may include determining optimum operating statuses (e.g., onor off) for a plurality of devices of the central plant equipment.According to an exemplary embodiment, the on/off combinations may bedetermined using binary optimization and quadratic compensation. Binaryoptimization may minimize a cost function representing the powerconsumption of devices in the applicable subplant. In some embodiments,non-exhaustive (i.e., not all potential combinations of devices areconsidered) binary optimization is used. Quadratic compensation may beused in considering devices whose power consumption is quadratic (andnot linear).

Step 1112 may include determining optimum operating setpoints forequipment using nonlinear optimization. Nonlinear optimization mayidentify operating setpoints that further minimize the low level costfunction J_(LL). In some embodiments, step 1112 includes providing theon/off decisions and setpoints to building automation system 108 for usein controlling the central plant equipment 60.

The construction and arrangement of the systems and methods as shown inthe various exemplary embodiments are illustrative only. Although only afew embodiments have been described in detail in this disclosure, manymodifications are possible (e.g., variations in sizes, dimensions,structures, shapes and proportions of the various elements, values ofparameters, mounting arrangements, use of materials, orientations,etc.). For example, the position of elements may be reversed orotherwise varied and the nature or number of discrete elements orpositions may be altered or varied. Accordingly, all such modificationsare intended to be included within the scope of the present disclosure.The order or sequence of any process or method steps may be varied orre-sequenced according to alternative embodiments. Other substitutions,modifications, changes, and omissions may be made in the design,operating conditions and arrangement of the exemplary embodimentswithout departing from the scope of the present disclosure.

The present disclosure contemplates methods, systems and programproducts on memory or other machine-readable media for accomplishingvarious operations. The embodiments of the present disclosure may beimplemented using existing computer processors, or by a special purposecomputer processor for an appropriate system, incorporated for this oranother purpose, or by a hardwired system. Embodiments within the scopeof the present disclosure include program products or memory comprisingmachine-readable media for carrying or having machine-executableinstructions or data structures stored thereon. Such machine-readablemedia can be any available media that can be accessed by a generalpurpose or special purpose computer or other machine with a processor.By way of example, such machine-readable media can comprise RAM, ROM,EPROM, EEPROM, CD-ROM or other optical disk storage, magnetic diskstorage or other magnetic storage devices, or any other medium which canbe used to carry or store desired program code in the form ofmachine-executable instructions or data structures and which can beaccessed by a general purpose or special purpose computer or othermachine with a processor. Combinations of the above are also includedwithin the scope of machine-readable media. Machine-executableinstructions include, for example, instructions and data which cause ageneral purpose computer, special purpose computer, or special purposeprocessing machines to perform a certain function or group of functions.

Although the figures may show a specific order of method steps, theorder of the steps may differ from what is depicted. Also two or moresteps may be performed concurrently or with partial concurrence. Suchvariation will depend on the software and hardware systems chosen and ondesigner choice. All such variations are within the scope of thedisclosure. Likewise, software implementations could be accomplishedwith standard programming techniques with rule based logic and otherlogic to accomplish the various connection steps, processing steps,comparison steps and decision steps.

What is claimed is:
 1. A controller for equipment that operate to provide heating or cooling to a building or campus, the controller comprising: a processing circuit configured to: obtain utility rate data indicating a price of one or more input resources consumed by the equipment at each of a plurality of time steps of an optimization period; obtain an objective function that expresses a total monetary cost of operating the equipment over the optimization period as a function of the utility rate data and an amount of the one or more input resources consumed by the equipment at each of the plurality of time steps; determine a relationship between (i) the amount of the one or more input resources consumed by the equipment and (ii) an amount of one or more output resources produced by the equipment as a result of consuming the one or more input resources; optimize the objective function over the optimization period, and subject to a constraint based on the relationship, to determine a distribution of production of the one or more output resources across the equipment at each of the plurality of time steps; and operate the equipment to achieve the distribution at each of the plurality of time steps.
 2. The controller of claim 1, wherein the processing circuit is configured to predict the demand for the one or more output resources using at least one of input from sensors within the building or campus or weather data for a location of the building or campus.
 3. The controller of claim 1, wherein the equipment comprise a plurality of subplants of a central plant and operate to provide the heating or cooling to the building or campus by heating or cooling a working fluid that is circulated to air handlers of the building or campus, wherein the one or more output resources comprise at least one of a cooled working fluid or a heated working fluid.
 4. The controller of claim 1, wherein the equipment comprise at least one of air handlers, dampers, or fans of the building or campus and operate to provide the heating or cooling to the building or campus by delivering heated or cooled air to one or more zones of the building or campus, wherein the one or more output resources comprise at least one of the heated air or the cooled air.
 5. The controller of claim 1, wherein: the relationship defines the amount of the one or more output resources produced by the equipment as a function of the amount of the one or more input resources consumed by the equipment; and the processing circuit is configured to use the relationship to predict the amount of the one or more resources consumed by the equipment to produce a corresponding amount of the one or more output resources.
 6. The controller of claim 1, wherein the processing circuit is configured to determine or update the relationship using experimental data obtained by operating the equipment.
 7. The controller of claim 1, wherein the objective function comprises cost variables representing a monetary cost associated with each of the one or more input resources consumed by the equipment at each of the plurality of time steps.
 8. A controller for equipment that operate to provide heating or cooling to a building or campus, the controller comprising: a processing circuit configured to: obtain utility rate data indicating a price of one or more input resources consumed by the equipment at each of a plurality of time steps of an optimization period; obtain an objective function that expresses a total monetary cost of operating the equipment over the optimization period as a function of the utility rate data and an amount of the one or more input resources consumed by the equipment at each of the plurality of time steps, wherein the objective function comprises: cost variables representing a monetary cost associated with each of the one or more input resources consumed by the equipment at each of the plurality of time steps; and production variables representing amounts of one or more output resources produced by the equipment as a result of consuming the one or more input resources, wherein the production variables are based on the amount of the one or more input resources consumed by the equipment at each of the plurality of time steps; optimize the objective function over the optimization period to determine a distribution of production of the one or more output resources across the equipment at each of the plurality of time steps, wherein determining the distribution comprises determining values for the production variables in the objective function; and operate the equipment to achieve the distribution at each of the plurality of time steps.
 9. The controller of claim 8, wherein the processing circuit is configured to predict demand for the one or more output resources using at least one of input from sensors within the building or campus or weather data for a location of the building or campus.
 10. The controller of claim 8, wherein the equipment comprise a plurality of subplants of a central plant and operate to provide the heating or cooling to the building or campus by heating or cooling a working fluid that is circulated to air handlers of the building or campus.
 11. The controller of claim 8, wherein the equipment comprise at least one of air handlers, dampers, or fans of the building or campus and operate to provide the heating or cooling to the building or campus by delivering heated or cooled air to one or more zones of the building or campus.
 12. The controller of claim 8, wherein the processing circuit is configured to: determine a relationship between the amount of the one or more input resources consumed by the equipment and the amount of the one or more output resources produced by the equipment; and optimize the objective function over the optimization period subject to a constraint, wherein the constraint is based on the relationship.
 13. A method for operating equipment to provide heating or cooling to a building or campus, the method comprising: obtaining utility rate data indicating a price of one or more input resources consumed by the equipment at each of a plurality of time steps of an optimization period; obtaining an objective function that expresses a total monetary cost of operating the equipment over the optimization period as a function of the utility rate data and an amount of the one or more input resources consumed by the equipment at each of the plurality of time steps; determining a relationship between (i) the amount of the one or more input resources consumed by the equipment and (ii) an amount of one or more output resources produced by the equipment as a result of consuming the one or more input resources; optimizing the objective function over the optimization period, and subject to a constraint based on the relationship, to determine a distribution of production of the one or more output resources across the equipment at each of the plurality of time steps; and operating the equipment to achieve the distribution at each of the plurality of time steps.
 14. The method of claim 13, comprising predicting the demand for the one or more output resources using at least one of input from sensors within the building or campus or weather data for a location of the building or campus.
 15. The method of claim 13, wherein the equipment comprise a plurality of subplants of a central plant and operate to provide the heating or cooling to the building or campus by heating or cooling a working fluid that is circulated to air handlers of the building or campus.
 16. The method of claim 13, wherein the equipment comprise at least one of air handlers, dampers, or fans of the building or campus and operate to provide the heating or cooling to the building or campus by delivering heated or cooled air to one or more zones of the building or campus.
 17. The method of claim 13, wherein the relationship defines the amount of the one or more output resources produced by the equipment as a function of the amount of the one or more input resources consumed by the equipment; the method comprising using the relationship to predict the amount of the one or more resources consumed by the equipment to produce a corresponding amount of the one or more output resources.
 18. The method of claim 13, comprising determining or updating the relationship using experimental data obtained by operating the equipment.
 19. The method of claim 13, wherein the objective function comprises: cost variables representing a monetary cost associated with each of the one or more input resources consumed by the equipment at each of the plurality of time steps; and production variables representing an energy load of the equipment at each of the plurality of time steps.
 20. The method of claim 13, comprising performing the optimization subject to a constraint that requires the amount of the one or more output resources produced by the equipment to satisfy demand of the building or campus for the one or more output resources. 